{
  "nbformat": 4,
  "nbformat_minor": 0,
  "metadata": {
    "kernelspec": {
      "display_name": "Python 3",
      "language": "python",
      "name": "python3"
    },
    "language_info": {
      "codemirror_mode": {
        "name": "ipython",
        "version": 3
      },
      "file_extension": ".py",
      "mimetype": "text/x-python",
      "name": "python",
      "nbconvert_exporter": "python",
      "pygments_lexer": "ipython3",
      "version": "3.7.6"
    },
    "colab": {
      "name": "Tool.ipynb",
      "provenance": [],
      "collapsed_sections": [
        "_iTuZpfIeXW1",
        "LgYLhVcDeXXA",
        "EXEZPdSfeXXC"
      ]
    }
  },
  "cells": [
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "OfsTg1f_-Pq4",
        "colab_type": "text"
      },
      "source": [
        "**Tasks Accomplished**\n",
        "\n",
        "1) Integrating my program with Kartikaya's more smoothly\n",
        "\n",
        "    a) Figuring out the use of the synthetic control models\n",
        "    b) Adapting the code to the various changes in variable names, format of data, etc.\n",
        "\n",
        "2) Predictor variables have been used to find the similar counties (this has yielded better results than just lagged-outcome variables).\n",
        "\n",
        "    a) The predictor variables used are population, median age, and population density. \n",
        "    Population and median age are drawn from the Census Bureau. \n",
        "    To get population density, I webscraped a wikipedia article and extracted \n",
        "    the land areas of the counties. I then divided the poplation column by the land area column.\n",
        "\n",
        "    b) This was done because the lagged-outcome-variable model yielded incorrect results. \n",
        "    The prediction not match the pre-intervention period and the post-intervention period \n",
        "    kept showing that masks were driving up cases.\n",
        "         i) This was because the donor pool counties had a much lower population (and presumably \n",
        "         population density), even though the initial number of cases matched my selected county.\n",
        "\n",
        "    c) With a combination of lagged-outcome variables and predictor variables,\n",
        "    the model is working well.\n",
        "\n",
        "**Remaining Tasks**\n",
        "\n",
        "1) Do multiple counties simultaneously.\n",
        "\n",
        "2) Look at the counties of states within the same region instead of limiting ourselves to the counties of any given state\n",
        "\n",
        "**Long-term**\n",
        "1) TEX, LATEX, and Overleaf"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "EAPN3HEZeXWE",
        "colab_type": "text"
      },
      "source": [
        "## Importing Libraries "
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "ot0eI208_We_",
        "colab_type": "code",
        "cellView": "both",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 34
        },
        "outputId": "38600ecf-dfe6-42b0-9bea-3ffbd0cc31c1"
      },
      "source": [
        "#@title DO NOT RUN THIS CELL UNLESS YOU HAVE BEEN DISCONNECTED FROM A SESSION\n",
        "from google.colab import drive\n",
        "drive.mount('/content/gdrive') \n",
        "#Mounts our drive onto our google colab system"
      ],
      "execution_count": 178,
      "outputs": [
        {
          "output_type": "stream",
          "text": [
            "Drive already mounted at /content/gdrive; to attempt to forcibly remount, call drive.mount(\"/content/gdrive\", force_remount=True).\n"
          ],
          "name": "stdout"
        }
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "rNnvqjTy4pT5",
        "colab_type": "code",
        "cellView": "both",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 71
        },
        "outputId": "165288e7-5ca6-4090-8a4f-a51a068fe8f9"
      },
      "source": [
        "#@title DO NOT RUN THIS CELL WITHOUT FIRST RESTARTING RUNTIME\n",
        "import sys\n",
        "!pip uninstall tslib\n",
        "#There are two modules called tslib. The one that !pip install accesses is NOT the one we want\n",
        "#This code will uninstall the irrelevant tslib library so that we can then access the relevant one\n",
        "\n",
        "sys.path.append('/content/gdrive/Shared drives/dev-SyntheticControl/COVID19-synthetic-control-analysis')\n",
        "print(sys.path)\n",
        "#Adds the synthetic control library from the google drive to the google colab system"
      ],
      "execution_count": 179,
      "outputs": [
        {
          "output_type": "stream",
          "text": [
            "\u001b[33mWARNING: Skipping tslib as it is not installed.\u001b[0m\n",
            "['', '/env/python', '/usr/lib/python36.zip', '/usr/lib/python3.6', '/usr/lib/python3.6/lib-dynload', '/usr/local/lib/python3.6/dist-packages', '/usr/lib/python3/dist-packages', '/usr/local/lib/python3.6/dist-packages/IPython/extensions', '/root/.ipython', '/content/gdrive/Shared drives/dev-SyntheticControl/COVID19-synthetic-control-analysis', '/content/gdrive/Shared drives/dev-SyntheticControl/COVID19-synthetic-control-analysis']\n"
          ],
          "name": "stdout"
        }
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "EfyONqnKtoh6",
        "colab_type": "code",
        "colab": {}
      },
      "source": [
        "#Synthetic control\n",
        "import tslib\n",
        "import tslib.src\n",
        "from tslib.src import tsUtils\n",
        "from tslib.tests import testdata\n",
        "from tslib.src.synthcontrol.syntheticControl import RobustSyntheticControl\n",
        "\n",
        "#Matplotlib\n",
        "import matplotlib.pyplot as plt\n",
        "import matplotlib.ticker as ticker\n",
        "import matplotlib.dates as mdates\n",
        "import matplotlib as mpl\n",
        "\n",
        "#Webscraping\n",
        "import requests\n",
        "from bs4 import BeautifulSoup\n",
        "\n",
        "#Miscellaneous\n",
        "import seaborn as sns\n",
        "import copy\n",
        "import pandas as pd\n",
        "import numpy as np\n",
        "from datetime import datetime as dtm\n",
        "from datetime import timedelta as td"
      ],
      "execution_count": 180,
      "outputs": []
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "iNAwCdkWeXWK",
        "colab_type": "code",
        "cellView": "form",
        "colab": {}
      },
      "source": [
        "#@title Copy-pasted Code from the Synthetic Control paper (not necessary to run)\n",
        "\n",
        "################################################################\n",
        "#\n",
        "# Robust Synthetic Control \n",
        "#\n",
        "# Implementation based on: \n",
        "# http://www.jmlr.org/papers/volume19/17-777/17-777.pdf\n",
        "#\n",
        "################################################################\n",
        "import numpy as np\n",
        "import pandas as pd\n",
        "\n",
        "from tslib.src.models.tsSVDModel import SVDModel\n",
        "from tslib.src.models.tsALSModel import ALSModel\n",
        "from tslib.src import tsUtils\n",
        "\n",
        "class RobustSyntheticControl(object):\n",
        "    \n",
        "    # seriesToPredictKey:       (string) the series of interest (key)\n",
        "    # kSingularValuesToKeep:    (int) the number of singular values to retain\n",
        "    # M:                        (int) the number of columns for the matrix\n",
        "    # probObservation:          (float) the independent probability of observation of each entry in the matrix\n",
        "    # modelType:                (string) SVD or ALS. Default is \"SVD\"\n",
        "    # svdMethod:                (string) the SVD method to use (optional)\n",
        "    # otherSeriesKeysArray:     (array) an array of keys for other series which will be used to predict \n",
        "\n",
        "    def __init__(self, seriesToPredictKey, kSingularValuesToKeep, M, probObservation=1.0, modelType='svd', svdMethod='numpy', otherSeriesKeysArray=[]):\n",
        "\n",
        "        self.seriesToPredictKey = seriesToPredictKey\n",
        "        self.otherSeriesKeysArray = otherSeriesKeysArray\n",
        "        self.N = 1 # each series is on its own row\n",
        "        self.M = M\n",
        "        self.kSingularValues = kSingularValuesToKeep\n",
        "        self.svdMethod = svdMethod\n",
        "        self.p = probObservation\n",
        "\n",
        "        if (modelType == 'als'):\n",
        "            self.model = ALSModel(self.seriesToPredictKey, self.kSingularValues, self.N, self.M, probObservation=self.p, otherSeriesKeysArray=self.otherSeriesKeysArray, includePastDataOnly=False)\n",
        "        elif (modelType == 'svd'):\n",
        "            self.model = SVDModel(self.seriesToPredictKey, self.kSingularValues, self.N, self.M, probObservation=self.p, svdMethod='numpy', otherSeriesKeysArray=self.otherSeriesKeysArray, includePastDataOnly=False)\n",
        "        else:\n",
        "            self.model = SVDModel(self.seriesToPredictKey, self.kSingularValues, self.N, self.M, probObservation=self.p, svdMethod='numpy', otherSeriesKeysArray=self.otherSeriesKeysArray, includePastDataOnly=False)\n",
        "\n",
        "        self.control = None # these are the synthetic control weights\n",
        "\n",
        "\n",
        "    # keyToSeriesDictionary: (Pandas dataframe) a key-value Series\n",
        "    # Note that the keys provided in the constructor MUST all be present\n",
        "    # The values must be all numpy arrays of floats.\n",
        "    def fit(self, keyToSeriesDF):\n",
        "        self.model.fit(keyToSeriesDF)\n",
        "\n",
        "\n",
        "    # otherKeysToSeriesDFNew:     (Pandas dataframe) needs to contain all keys provided in the model;\n",
        "    #                               all series/array MUST be of length >= 1, \n",
        "    #                               If longer than 1, then the most recent point will be used (for each series)\n",
        "    def predict(self, otherKeysToSeriesDFNew):\n",
        "        prediction = np.dot(self.model.weights, otherKeysToSeriesDFNew[self.otherSeriesKeysArray].T)\n",
        "        return prediction\n",
        "\n",
        "    # return the synthetic control weights\n",
        "    def getControl():\n",
        "        if (self.model.weights is None):\n",
        "            raise Exception('Before calling getControl() you need to call the fit() method first.')\n",
        "        else:\n",
        "            return self.model.weights"
      ],
      "execution_count": 181,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "Hifrik7neXWM",
        "colab_type": "text"
      },
      "source": [
        "## Input Variables "
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "GZ11Y8uJeXWN",
        "colab_type": "code",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 52
        },
        "outputId": "dd271e3d-b81d-40c6-a201-d7c89cbba141"
      },
      "source": [
        "'''\n",
        "format = \"%Y-%m-%d\"\n",
        "\n",
        "target=input('Enter the state: ')\n",
        "county=input('Enter the county in that state: ')\n",
        "date=input('Enter the intervention date in '+format+' format: ')\n",
        "start=input('Starting date in '+format+' format: ')\n",
        "end=input('Ending date in '+format+' format: ')\n",
        "'''"
      ],
      "execution_count": 182,
      "outputs": [
        {
          "output_type": "execute_result",
          "data": {
            "application/vnd.google.colaboratory.intrinsic+json": {
              "type": "string"
            },
            "text/plain": [
              "'\\nformat = \"%Y-%m-%d\"\\n\\ntarget=input(\\'Enter the state: \\')\\ncounty=input(\\'Enter the county in that state: \\')\\ndate=input(\\'Enter the intervention date in \\'+format+\\' format: \\')\\nstart=input(\\'Starting date in \\'+format+\\' format: \\')\\nend=input(\\'Ending date in \\'+format+\\' format: \\')\\n'"
            ]
          },
          "metadata": {
            "tags": []
          },
          "execution_count": 182
        }
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "9tOyDmfVeXWW",
        "colab_type": "code",
        "colab": {}
      },
      "source": [
        "URL1 ='https://raw.githubusercontent.com/nytimes/covid-19-data/master/us-counties.csv'\n",
        "URL2 = \"https://raw.githubusercontent.com/nytimes/covid-19-data/master/mask-use/mask-use-by-county.csv\"\n",
        "URL3 = \"https://www.census.gov/geographies/reference-files/2017/demo/popest/2017-fips.html\"\n",
        "URL4 = \"https://raw.githubusercontent.com/jasonong/List-of-US-States/master/states.csv\""
      ],
      "execution_count": 183,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "43Q53bBEeXWf",
        "colab_type": "text"
      },
      "source": [
        "## Finding Similar Counties "
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "yS62IZ2Mgtr7",
        "colab_type": "code",
        "cellView": "form",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 1000
        },
        "outputId": "a42351bf-517f-4665-b2d2-3edc1f67f9f7"
      },
      "source": [
        "# @title No need to run this for now\n",
        "nyTimesData = pd.read_csv(\"https://raw.githubusercontent.com/nytimes/covid-19-data/master/us-counties.csv\") #gets the county data\n",
        "nyTimesDataAdded2 = pd.DataFrame()\n",
        "\n",
        "for state_name, state in nyTimesData.groupby(\"state\"): #iterates through each state\n",
        "  nyTimesDataAdded = pd.DataFrame()\n",
        "  for county_name, county in state.groupby(\"county\"): #iterates through each county of each state\n",
        "    county[\"new cases\"] = county[\"cases\"].diff() #gets the new cases per day\n",
        "    county[\"new deaths\"] = county[\"deaths\"].diff() #gets the new deaths per day\n",
        "    nyTimesDataAdded = pd.concat([nyTimesDataAdded, county]) #adds the information about the county to the overall dataframe about the state\n",
        "  \n",
        "  nyTimesDataAdded2 = pd.concat([nyTimesDataAdded2, nyTimesDataAdded]) #adds the information about the state to the overall dataframe about the county\n",
        "  print(state_name)\n",
        "\n",
        "nyTimesData = nyTimesDataAdded2.reset_index(drop = True).set_index([\"date\", \"county\", \"state\"]).drop(\"fips\", axis = 1)\n",
        "nyTimesData"
      ],
      "execution_count": 184,
      "outputs": [
        {
          "output_type": "stream",
          "text": [
            "/usr/local/lib/python3.6/dist-packages/ipykernel_launcher.py:8: SettingWithCopyWarning: \n",
            "A value is trying to be set on a copy of a slice from a DataFrame.\n",
            "Try using .loc[row_indexer,col_indexer] = value instead\n",
            "\n",
            "See the caveats in the documentation: https://pandas.pydata.org/pandas-docs/stable/user_guide/indexing.html#returning-a-view-versus-a-copy\n",
            "  \n",
            "/usr/local/lib/python3.6/dist-packages/ipykernel_launcher.py:9: SettingWithCopyWarning: \n",
            "A value is trying to be set on a copy of a slice from a DataFrame.\n",
            "Try using .loc[row_indexer,col_indexer] = value instead\n",
            "\n",
            "See the caveats in the documentation: https://pandas.pydata.org/pandas-docs/stable/user_guide/indexing.html#returning-a-view-versus-a-copy\n",
            "  if __name__ == '__main__':\n"
          ],
          "name": "stderr"
        },
        {
          "output_type": "stream",
          "text": [
            "Alabama\n",
            "Alaska\n",
            "Arizona\n",
            "Arkansas\n",
            "California\n",
            "Colorado\n",
            "Connecticut\n",
            "Delaware\n",
            "District of Columbia\n",
            "Florida\n",
            "Georgia\n",
            "Guam\n",
            "Hawaii\n",
            "Idaho\n",
            "Illinois\n",
            "Indiana\n",
            "Iowa\n",
            "Kansas\n",
            "Kentucky\n",
            "Louisiana\n",
            "Maine\n",
            "Maryland\n",
            "Massachusetts\n",
            "Michigan\n",
            "Minnesota\n",
            "Mississippi\n",
            "Missouri\n",
            "Montana\n",
            "Nebraska\n",
            "Nevada\n",
            "New Hampshire\n",
            "New Jersey\n",
            "New Mexico\n",
            "New York\n",
            "North Carolina\n",
            "North Dakota\n",
            "Northern Mariana Islands\n",
            "Ohio\n",
            "Oklahoma\n",
            "Oregon\n",
            "Pennsylvania\n",
            "Puerto Rico\n",
            "Rhode Island\n",
            "South Carolina\n",
            "South Dakota\n",
            "Tennessee\n",
            "Texas\n",
            "Utah\n",
            "Vermont\n",
            "Virgin Islands\n",
            "Virginia\n",
            "Washington\n",
            "West Virginia\n",
            "Wisconsin\n",
            "Wyoming\n"
          ],
          "name": "stdout"
        },
        {
          "output_type": "execute_result",
          "data": {
            "text/html": [
              "<div>\n",
              "<style scoped>\n",
              "    .dataframe tbody tr th:only-of-type {\n",
              "        vertical-align: middle;\n",
              "    }\n",
              "\n",
              "    .dataframe tbody tr th {\n",
              "        vertical-align: top;\n",
              "    }\n",
              "\n",
              "    .dataframe thead th {\n",
              "        text-align: right;\n",
              "    }\n",
              "</style>\n",
              "<table border=\"1\" class=\"dataframe\">\n",
              "  <thead>\n",
              "    <tr style=\"text-align: right;\">\n",
              "      <th></th>\n",
              "      <th></th>\n",
              "      <th></th>\n",
              "      <th>cases</th>\n",
              "      <th>deaths</th>\n",
              "      <th>new cases</th>\n",
              "      <th>new deaths</th>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>date</th>\n",
              "      <th>county</th>\n",
              "      <th>state</th>\n",
              "      <th></th>\n",
              "      <th></th>\n",
              "      <th></th>\n",
              "      <th></th>\n",
              "    </tr>\n",
              "  </thead>\n",
              "  <tbody>\n",
              "    <tr>\n",
              "      <th>2020-03-24</th>\n",
              "      <th>Autauga</th>\n",
              "      <th>Alabama</th>\n",
              "      <td>1</td>\n",
              "      <td>0</td>\n",
              "      <td>NaN</td>\n",
              "      <td>NaN</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>2020-03-25</th>\n",
              "      <th>Autauga</th>\n",
              "      <th>Alabama</th>\n",
              "      <td>4</td>\n",
              "      <td>0</td>\n",
              "      <td>3.0</td>\n",
              "      <td>0.0</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>2020-03-26</th>\n",
              "      <th>Autauga</th>\n",
              "      <th>Alabama</th>\n",
              "      <td>6</td>\n",
              "      <td>0</td>\n",
              "      <td>2.0</td>\n",
              "      <td>0.0</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>2020-03-27</th>\n",
              "      <th>Autauga</th>\n",
              "      <th>Alabama</th>\n",
              "      <td>6</td>\n",
              "      <td>0</td>\n",
              "      <td>0.0</td>\n",
              "      <td>0.0</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>2020-03-28</th>\n",
              "      <th>Autauga</th>\n",
              "      <th>Alabama</th>\n",
              "      <td>6</td>\n",
              "      <td>0</td>\n",
              "      <td>0.0</td>\n",
              "      <td>0.0</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>...</th>\n",
              "      <th>...</th>\n",
              "      <th>...</th>\n",
              "      <td>...</td>\n",
              "      <td>...</td>\n",
              "      <td>...</td>\n",
              "      <td>...</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>2020-09-03</th>\n",
              "      <th>Weston</th>\n",
              "      <th>Wyoming</th>\n",
              "      <td>20</td>\n",
              "      <td>0</td>\n",
              "      <td>1.0</td>\n",
              "      <td>0.0</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>2020-09-04</th>\n",
              "      <th>Weston</th>\n",
              "      <th>Wyoming</th>\n",
              "      <td>20</td>\n",
              "      <td>0</td>\n",
              "      <td>0.0</td>\n",
              "      <td>0.0</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>2020-09-05</th>\n",
              "      <th>Weston</th>\n",
              "      <th>Wyoming</th>\n",
              "      <td>20</td>\n",
              "      <td>0</td>\n",
              "      <td>0.0</td>\n",
              "      <td>0.0</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>2020-09-06</th>\n",
              "      <th>Weston</th>\n",
              "      <th>Wyoming</th>\n",
              "      <td>20</td>\n",
              "      <td>0</td>\n",
              "      <td>0.0</td>\n",
              "      <td>0.0</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>2020-09-07</th>\n",
              "      <th>Weston</th>\n",
              "      <th>Wyoming</th>\n",
              "      <td>20</td>\n",
              "      <td>0</td>\n",
              "      <td>0.0</td>\n",
              "      <td>0.0</td>\n",
              "    </tr>\n",
              "  </tbody>\n",
              "</table>\n",
              "<p>511829 rows × 4 columns</p>\n",
              "</div>"
            ],
            "text/plain": [
              "                            cases  deaths  new cases  new deaths\n",
              "date       county  state                                        \n",
              "2020-03-24 Autauga Alabama      1       0        NaN         NaN\n",
              "2020-03-25 Autauga Alabama      4       0        3.0         0.0\n",
              "2020-03-26 Autauga Alabama      6       0        2.0         0.0\n",
              "2020-03-27 Autauga Alabama      6       0        0.0         0.0\n",
              "2020-03-28 Autauga Alabama      6       0        0.0         0.0\n",
              "...                           ...     ...        ...         ...\n",
              "2020-09-03 Weston  Wyoming     20       0        1.0         0.0\n",
              "2020-09-04 Weston  Wyoming     20       0        0.0         0.0\n",
              "2020-09-05 Weston  Wyoming     20       0        0.0         0.0\n",
              "2020-09-06 Weston  Wyoming     20       0        0.0         0.0\n",
              "2020-09-07 Weston  Wyoming     20       0        0.0         0.0\n",
              "\n",
              "[511829 rows x 4 columns]"
            ]
          },
          "metadata": {
            "tags": []
          },
          "execution_count": 184
        }
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "_j-Tb5VssPEd",
        "colab_type": "code",
        "colab": {}
      },
      "source": [
        "def get_raw_data(URL1, ST_NAME):\n",
        "  state = pd.read_csv(URL1, index_col=0)\n",
        "  looking_for_alternative = \"cases\"\n",
        "  state=state[state['state']==target].pivot_table(index='county', columns='date', values= looking_for_alternative).fillna(0)\n",
        "  state_data_of_interest = (pd.read_csv(URL1).\n",
        "                            set_index([\"date\", \"county\", \"state\"]).\n",
        "                            drop(\"fips\", axis = 1).\n",
        "                            xs(ST_NAME, level = \"state\").\n",
        "                            reset_index().\n",
        "                            sort_values(by = [\"date\", \"county\"]).\n",
        "                            set_index([\"date\", \"county\"]))\n",
        "\n",
        "  dates = []\n",
        "  column_names = []\n",
        "  for i in days_prior:\n",
        "    dates.append(dtm.strftime(dtm.strptime(intervention, format) - td(days = i), format))\n",
        "    if i == 1:\n",
        "      column_names.append(\"Scaled magnitude of \"+looking_for_alternative+\" (\"+str(i)+\" day prior)\")\n",
        "    else:\n",
        "      column_names.append(\"Scaled magnitude of \"+looking_for_alternative+\" (\"+str(i)+\" days prior)\")\n",
        "  raw_X = (state[dates].\n",
        "      rename(columns = dict(zip(state[dates].columns, \n",
        "                                column_names))).\n",
        "      replace({0:np.nan}).\n",
        "      dropna())\n",
        "  return raw_X, raw_X.applymap(np.log)"
      ],
      "execution_count": 426,
      "outputs": []
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "nhWuIFvClU4m",
        "colab_type": "code",
        "colab": {}
      },
      "source": [
        "def get_age_database_links(URL4):\n",
        "\n",
        "  URL4 = \"https://www.census.gov/data/tables/time-series/demo/popest/2010s-counties-detail.html\" #link for webscraping\n",
        "  a = BeautifulSoup(requests.get(URL4).text, \"html.parser\")\n",
        "\n",
        "  age_html = [i for i in [i for i in  [i for i in [i for i in a.find_all(\"div\") if len(i.get_attribute_list(\"class\"))>0 ]]  if not( i.get_attribute_list(\"class\")[0] in [None, \"\"])] if i.get_attribute_list(\"class\")[0].startswith(\"state\")][0]\n",
        "  age_databases = pd.DataFrame( data = [[i.get_text().strip(), i.get_attribute_list(\"href\")[0] ] for i in age_html.find_all(\"a\")], columns = [\"State\", \"Link\"]).set_index(\"State\")\n",
        "  #webscrapes the URL and gets the links to the csv files that contain the census information for each state\n",
        "\n",
        "  return age_databases"
      ],
      "execution_count": 186,
      "outputs": []
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "QbxZso4DmHsw",
        "colab_type": "code",
        "colab": {}
      },
      "source": [
        "def add_population_to_X(X, age_databases):\n",
        "\n",
        "  main_df = pd.DataFrame()\n",
        "  for state in age_databases.index: #iterates through the states\n",
        "    if state == target:\n",
        "      try:\n",
        "        x = pd.read_csv(\"https:\" + age_databases[\"Link\"].loc[state]).drop([\"SUMLEV\", \"STNAME\", \"COUNTY\"], axis = 1)\n",
        "        #reads the csv file for each state\n",
        "\n",
        "        x = x.groupby(\"CTYNAME\").apply(np.mean).drop([\"YEAR\"], axis = 1).reset_index()\n",
        "        #the csv file has information for six years about each county\n",
        "        #this averages together the census data across those six years for each county\n",
        "        #the dataframe is grouped by county and the mean function is applied to it\n",
        "\n",
        "        x[\"STATE\"] = state\n",
        "\n",
        "        land_area = pd.read_html(\"https://en.wikipedia.org/wiki/List_of_counties_in_\"+state.replace(\" \", \"_\"))\n",
        "        #for webscraping the area of each county in the state (allowing me to calculate the population density of each county)\n",
        "        land_area = [i for i in land_area if len(i.columns) > 8][0]\n",
        "        #goes to the correct table on the wikipedia page\n",
        "\n",
        "        division = land_area.columns[0] #this is the county/parish\n",
        "        land_area = land_area.drop([\"Map\"], axis = 1, errors = \"ignore\").set_index(division)\n",
        "        #Most wikipedia pages have a column that shows a map of where the county is; I'm dropping it\n",
        "        #I am ignoring errors, because those will probably come from pages that don't have that column, in which case nothing needs to be done\n",
        "        #the county is set as the index\n",
        "\n",
        "        land_area = land_area[ [ land_area.columns[-1]] ].applymap(lambda area: float(\"\".join((area.split()[0]).split(\",\")))  ).rename(columns = {land_area.columns[-1]:\"Land Area\"}).applymap(round).reset_index()\n",
        "        #fixes the formatting of the numbers in the area column (they are strings with commas ever three digits)\n",
        "        \n",
        "        x = x.merge(land_area, left_on = \"CTYNAME\", right_on = division).drop([division], axis = 1)\n",
        "        #merges the census dataframe and the wikipedia dataframe into a single dataframe about the state\n",
        "\n",
        "        x[\"Population Density\"] = x[\"POPESTIMATE\"] / x[\"Land Area\"]\n",
        "        main_df = pd.concat([main_df, x])\n",
        "        #adds in the dataframe about the state to the overall dataframe about the country\n",
        "        print(state)\n",
        "      except:\n",
        "        print(state + \" (Could not add)\") #If anything went wrong, I'm not sure how to fix it, so I'll just ignore that state\n",
        "  main_df = (main_df.rename(columns = {\"CTYNAME\":\"County\", \"STATE\":\"State\"}).\n",
        "            set_index([\"State\", \"County\"])\n",
        "            [[\"POPESTIMATE\", \"MEDIAN_AGE_TOT\", \"Population Density\"]].\n",
        "            rename(columns = dict(zip ( [\"POPESTIMATE\", \"MEDIAN_AGE_TOT\"],\n",
        "                                        [\"Population\", \"Median Age\"]) ),\n",
        "                    errors = \"ignore\" ).\n",
        "            xs(target))\n",
        "\n",
        "  main_df.index = [i.split(\" County\")[0] for i in main_df.index]\n",
        "  X = X.join(main_df)\n",
        "  return X"
      ],
      "execution_count": 187,
      "outputs": []
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "d9KEKJDQsc4H",
        "colab_type": "code",
        "colab": {}
      },
      "source": [
        "def further_scaling(X):\n",
        "  X[\"Population\"] = X[\"Population\"].apply(np.log)\n",
        "  X[\"Population Density\"] = X[\"Population Density\"].apply(np.log)\n",
        "\n",
        "  for i in X.columns:\n",
        "    X[i] = 100*(X[i] - min(X[i]) ) / (max(X[i]) - min(X[i])) #min-max scaling\n",
        "\n",
        "  X.columns = list(X.columns[:-3])+[i+\" (Scaled)\" for i in X.columns[-3:] ]\n",
        "  return X"
      ],
      "execution_count": 188,
      "outputs": []
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "ZZJqBeo3gwY1",
        "colab_type": "code",
        "cellView": "form",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 95
        },
        "outputId": "f03d6b9f-703e-4e85-cac8-b192d8694157"
      },
      "source": [
        "#@title No need to run this now\n",
        "for i in days_prior: #Gets the Covid-19 information from some days before the intervention date\n",
        "  date = dtm.strftime(dtm.strptime(intervention, format) - td(days = i), format) #Gets the date i days before the intervention date\n",
        "  state_data_of_interest_subset = state_data_of_interest.xs(date, level = \"date\").reset_index() #Gets the information from that date\n",
        "  state_data_of_interest_subset[\"county\"] += \" County\" #Changes the name of the county to match the other dataframe (e.g. Allen --> Allen County; Shawnee --> Shawnee County)\n",
        "  state_data_of_interest_subset = state_data_of_interest_subset.set_index(\"county\") #Sets the county as the index (helps to join the dataframes)\n",
        "  if i==1:\n",
        "    state_data_of_interest_subset.columns = [j + \" (\"+str(i)+\" day prior)\" for j in state_data_of_interest_subset.columns]\n",
        "  else:\n",
        "    state_data_of_interest_subset.columns = [j + \" (\"+str(i)+\" days prior)\" for j in state_data_of_interest_subset.columns]\n",
        "  \n",
        "  state_data_of_interest_subset = state_data_of_interest_subset[state_data_of_interest_subset.columns[:-1]] #gets rid of deaths per day\n",
        "  #I did this because the number of deaths per day is always either 0 or 1 (with a few exceptions)\n",
        "  #This minimal difference gets amplified by the min-max scaling to the point where I feel it has too much influence\n",
        "  #Without the min-max scaling, it has too little influence to be relevant\n",
        "\n",
        "  X = X.join(state_data_of_interest_subset, how = \"inner\") #joins the dataframes\n",
        "\n",
        "X = X.fillna(0) #The number of new cases and new deaths on the first day is NaN; this will replace these values with zero\n",
        "\n",
        "\n",
        "X = X[[i for i in X.columns if (\"cases\" in i) and not(\"new\" in i)]]\n",
        "Y = X\n",
        "X = X.applymap(np.log)\n",
        "\n",
        "\n",
        "#for i in X.columns:\n",
        "#  X[i] = 100*(X[i] - min(X[i]) ) / (max(X[i]) - min(X[i])) #min-max scaling\n",
        "X"
      ],
      "execution_count": 189,
      "outputs": [
        {
          "output_type": "execute_result",
          "data": {
            "text/html": [
              "<div>\n",
              "<style scoped>\n",
              "    .dataframe tbody tr th:only-of-type {\n",
              "        vertical-align: middle;\n",
              "    }\n",
              "\n",
              "    .dataframe tbody tr th {\n",
              "        vertical-align: top;\n",
              "    }\n",
              "\n",
              "    .dataframe thead th {\n",
              "        text-align: right;\n",
              "    }\n",
              "</style>\n",
              "<table border=\"1\" class=\"dataframe\">\n",
              "  <thead>\n",
              "    <tr style=\"text-align: right;\">\n",
              "      <th></th>\n",
              "      <th>Scaled magnitude of cases (1 day prior)</th>\n",
              "      <th>Scaled magnitude of cases (7 days prior)</th>\n",
              "      <th>Scaled magnitude of cases (15 days prior)</th>\n",
              "      <th>cases (1 day prior)</th>\n",
              "      <th>cases (7 days prior)</th>\n",
              "      <th>cases (15 days prior)</th>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>county</th>\n",
              "      <th></th>\n",
              "      <th></th>\n",
              "      <th></th>\n",
              "      <th></th>\n",
              "      <th></th>\n",
              "      <th></th>\n",
              "    </tr>\n",
              "  </thead>\n",
              "  <tbody>\n",
              "  </tbody>\n",
              "</table>\n",
              "</div>"
            ],
            "text/plain": [
              "Empty DataFrame\n",
              "Columns: [Scaled magnitude of cases (1 day prior), Scaled magnitude of cases (7 days prior), Scaled magnitude of cases (15 days prior), cases (1 day prior), cases (7 days prior), cases (15 days prior)]\n",
              "Index: []"
            ]
          },
          "metadata": {
            "tags": []
          },
          "execution_count": 189
        }
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "vzRbkTx5hUz2",
        "colab_type": "code",
        "colab": {}
      },
      "source": [
        "#Gets the New York Times data\n",
        "\n",
        "def NYTimesMaskData(URL2):\n",
        "\n",
        "  nyTimesMaskPolicy = pd.read_csv(URL2, dtype=str).set_index(\"COUNTYFP\")\n",
        "  #************************************************************************************************************************\n",
        "  a = dtm.now()\n",
        "  soup = BeautifulSoup(requests.get(URL3).text, \"html.parser\")\n",
        "  spreadsheets = [i for i in soup.find_all(\"a\") if i.get_attribute_list(\"class\")[0] == \"uscb-text-link\"]\n",
        "  countyFIPS_link = [i for i in spreadsheets if \"County\" in i.get_text()][0].get_attribute_list(\"href\")[0]\n",
        "  countyFIPS = pd.read_excel(\"https:\"+countyFIPS_link, skiprows = 4, dtype=str)\n",
        "  countyFIPS[\"County Code\"] = countyFIPS[\"State Code (FIPS)\"] + countyFIPS[\"County Code (FIPS)\"]\n",
        "  countyFIPS = countyFIPS.rename(columns = {\"Area Name (including legal/statistical area description)\":\"Place\"})[[\"Place\", \"County Code\", \"State Code (FIPS)\", \"County Code (FIPS)\"]]\n",
        "\n",
        "  states = list(pd.read_csv(\"https://raw.githubusercontent.com/jasonong/List-of-US-States/master/states.csv\")[\"State\"])\n",
        "  stateFIPS = countyFIPS[countyFIPS[\"Place\"].isin(states)].set_index(\"Place\").applymap(lambda x: x[:2])\n",
        "  countyFIPS = countyFIPS[countyFIPS[\"Place\"].str.lower().apply(lambda x: x.endswith(\"county\"))]\n",
        "  countyFIPS[\"State\"] = \"\"\n",
        "\n",
        "  for i in countyFIPS.index:\n",
        "    number = countyFIPS.loc[i, \"State Code (FIPS)\"]\n",
        "    countyFIPS.loc[i, \"State\"] = stateFIPS.reset_index().set_index(\"County Code\").loc[number][\"Place\"]\n",
        "\n",
        "  state_county_FIPS = (countyFIPS.query(\"State == '\"+ST_NAME+\"'\").\n",
        "                      reset_index(drop = True).\n",
        "                      applymap(lambda x: x.split(\" County\")[0]))\n",
        "\n",
        "  b = dtm.now()\n",
        "  print(\"Time to run function:\", b-a)\n",
        "  state_county_FIPS\n",
        "  #************************************************************************************************************************\n",
        "\n",
        "  weighted_sum = 0\n",
        "  for i in range(5):\n",
        "    weighted_sum += i*nyTimesMaskPolicy[nyTimesMaskPolicy.columns[i]].apply(float)\n",
        "\n",
        "\n",
        "  relevant_FIPS_codes = list((state_county_FIPS[\"State Code (FIPS)\"] + state_county_FIPS[\"County Code (FIPS)\"]).values)\n",
        "  state_counties_mask_usage = (weighted_sum.loc[relevant_FIPS_codes]\n",
        "                              .reset_index().set_index(\"COUNTYFP\")\n",
        "                              .join(state_county_FIPS.set_index(\"County Code\"))[[\"NEVER\", \"Place\"]]\n",
        "                              .rename(columns = {\"NEVER\":\"Usage Index\", \"PLACE\":\"County\"}))\n",
        "  \n",
        "  return state_counties_mask_usage.sort_values(by = \"Usage Index\").iloc[[0, -1]], state_counties_mask_usage"
      ],
      "execution_count": 190,
      "outputs": []
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "hec5IpppTClh",
        "colab_type": "code",
        "colab": {}
      },
      "source": [
        "def find_closest_counties(df, county):\n",
        "  new_df = pd.DataFrame(columns = [\"Euclidean Distance from \"+county])\n",
        "  for value in (df.index):\n",
        "    new_df.loc[value] = (((df.loc[county] - df.loc[value])**2).sum())**0.5\n",
        "  return new_df"
      ],
      "execution_count": 191,
      "outputs": []
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "Naexb35Qrvfj",
        "colab_type": "code",
        "colab": {}
      },
      "source": [
        "divisions = 3\n",
        "middle_counterfactual = 2\n",
        "\n",
        "def extract_donor_pool(closeness_scores2, state_counties_mask_usage, counties_in_donor_pool,\n",
        "                       divisions = 2, middle_counterfactual = 2):\n",
        "  initial_county_number = round(len(closeness_scores2)/2)\n",
        "  closeness_scores = (closeness_scores2.\n",
        "                      join(state_counties_mask_usage.reset_index().\n",
        "                          set_index(\"Place\")).\n",
        "                      dropna())\n",
        "\n",
        "\n",
        "  mean_score = np.mean(closeness_scores[\"Usage Index\"])\n",
        "  standard_deviation = np.std(closeness_scores[\"Usage Index\"])\n",
        "\n",
        "  if divisions == 3:\n",
        "    closeness_scores[\"Group\"] = closeness_scores[\"Usage Index\"] - mean_score\n",
        "    for i in closeness_scores.index:\n",
        "      if abs(closeness_scores.loc[i][\"Group\"])<=standard_deviation:\n",
        "        closeness_scores.loc[i, \"Group\"] = 1\n",
        "      else:\n",
        "        closeness_scores.loc[i, \"Group\"] = 2*int(closeness_scores.loc[i][\"Usage Index\"] > mean_score)\n",
        "  else:\n",
        "    closeness_scores[\"Group\"] = 2*((closeness_scores[\"Usage Index\"] >= mean_score).apply(int)) \n",
        "  closeness_scores\n",
        "\n",
        "  county_of_interest_group = closeness_scores.loc[county_of_interest][\"Group\"]\n",
        "  counterfactual_group = 2*(county_of_interest_group < 1)\n",
        "  if county_of_interest_group == 1:\n",
        "    counterfactual_group = middle_counterfactual\n",
        "\n",
        "  donor_pool_table = (closeness_scores.\n",
        "                      query(\"Group == \"+str(counterfactual_group)).\n",
        "                      sort_values(by = closeness_scores.columns[0]).\n",
        "                      head(counties_in_donor_pool))\n",
        "  donor_pool = list(donor_pool_table.index)\n",
        "  return donor_pool, closeness_scores, mean_score, standard_deviation"
      ],
      "execution_count": 332,
      "outputs": []
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "Yw9F2c-vQvin",
        "colab_type": "code",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 920
        },
        "outputId": "d5490e0a-8f0a-4ba4-c6fa-0fea8b8fa31b"
      },
      "source": [
        "format = \"%Y-%m-%d\"\n",
        "target = \"South Carolina\"\n",
        "county = \"Williamsburg\"\n",
        "\n",
        "#Dutchess County, New York\n",
        "#Sumter, South Carolina\n",
        "\n",
        "date   = \"2020-06-02\"\n",
        "start  = \"2020-03-10\"\n",
        "end    = \"2020-07-10\"\n",
        "\n",
        "ST_NAME = target\n",
        "intervention = date\n",
        "county_of_interest = str(county)\n",
        "#days_prior = [1, 3, 5, 7, 9, 11, 13, 15] #How many days before the intervention date do we want to consider?\n",
        "days_prior = [1, 7, 15] #How many days before the intervention date do we want to consider?\n",
        "max_euclidean_distance = 50\n",
        "counties_in_donor_pool = 4\n",
        "raw_X, X = get_raw_data(URL1, ST_NAME)\n",
        "age_databases = get_age_database_links(URL4)\n",
        "X = add_population_to_X(X, age_databases)\n",
        "X = further_scaling(X)\n",
        "extremes, state_counties_mask_usage = NYTimesMaskData(URL2)\n",
        "\n",
        "closeness_scores2 = (find_closest_counties(X, county_of_interest).\n",
        "                    sort_values(\"Euclidean Distance from \"+county_of_interest).\n",
        "                    query(\"`Euclidean Distance from \"+county_of_interest+\"` <= \"+str(max_euclidean_distance)))\n",
        "\n",
        "donor_pool, closeness_scores, mean_score, standard_deviation = extract_donor_pool(closeness_scores2, state_counties_mask_usage, counties_in_donor_pool,\n",
        "                                                                                  divisions = 2, middle_counterfactual = 0)\n",
        "\n",
        "\n",
        "closeness_scores"
      ],
      "execution_count": 413,
      "outputs": [
        {
          "output_type": "stream",
          "text": [
            "South Carolina\n",
            "Time to run function: 0:00:08.962903\n"
          ],
          "name": "stdout"
        },
        {
          "output_type": "execute_result",
          "data": {
            "text/html": [
              "<div>\n",
              "<style scoped>\n",
              "    .dataframe tbody tr th:only-of-type {\n",
              "        vertical-align: middle;\n",
              "    }\n",
              "\n",
              "    .dataframe tbody tr th {\n",
              "        vertical-align: top;\n",
              "    }\n",
              "\n",
              "    .dataframe thead th {\n",
              "        text-align: right;\n",
              "    }\n",
              "</style>\n",
              "<table border=\"1\" class=\"dataframe\">\n",
              "  <thead>\n",
              "    <tr style=\"text-align: right;\">\n",
              "      <th></th>\n",
              "      <th>Euclidean Distance from Williamsburg</th>\n",
              "      <th>COUNTYFP</th>\n",
              "      <th>Usage Index</th>\n",
              "      <th>Group</th>\n",
              "    </tr>\n",
              "  </thead>\n",
              "  <tbody>\n",
              "    <tr>\n",
              "      <th>Williamsburg</th>\n",
              "      <td>0.000000</td>\n",
              "      <td>45089</td>\n",
              "      <td>3.295</td>\n",
              "      <td>2</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>Colleton</th>\n",
              "      <td>16.709852</td>\n",
              "      <td>45029</td>\n",
              "      <td>3.129</td>\n",
              "      <td>2</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>Saluda</th>\n",
              "      <td>17.057280</td>\n",
              "      <td>45081</td>\n",
              "      <td>3.063</td>\n",
              "      <td>2</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>Clarendon</th>\n",
              "      <td>17.925122</td>\n",
              "      <td>45027</td>\n",
              "      <td>2.818</td>\n",
              "      <td>0</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>Lee</th>\n",
              "      <td>18.898024</td>\n",
              "      <td>45061</td>\n",
              "      <td>2.863</td>\n",
              "      <td>0</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>Fairfield</th>\n",
              "      <td>18.935085</td>\n",
              "      <td>45039</td>\n",
              "      <td>3.318</td>\n",
              "      <td>2</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>Chesterfield</th>\n",
              "      <td>20.095427</td>\n",
              "      <td>45025</td>\n",
              "      <td>2.866</td>\n",
              "      <td>0</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>Marlboro</th>\n",
              "      <td>21.920866</td>\n",
              "      <td>45069</td>\n",
              "      <td>2.840</td>\n",
              "      <td>0</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>Chester</th>\n",
              "      <td>23.574235</td>\n",
              "      <td>45023</td>\n",
              "      <td>2.865</td>\n",
              "      <td>0</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>Edgefield</th>\n",
              "      <td>24.201982</td>\n",
              "      <td>45037</td>\n",
              "      <td>3.021</td>\n",
              "      <td>2</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>Barnwell</th>\n",
              "      <td>24.931590</td>\n",
              "      <td>45011</td>\n",
              "      <td>2.768</td>\n",
              "      <td>0</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>Abbeville</th>\n",
              "      <td>25.399076</td>\n",
              "      <td>45001</td>\n",
              "      <td>2.570</td>\n",
              "      <td>0</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>Newberry</th>\n",
              "      <td>25.962924</td>\n",
              "      <td>45071</td>\n",
              "      <td>2.974</td>\n",
              "      <td>2</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>Marion</th>\n",
              "      <td>26.588717</td>\n",
              "      <td>45067</td>\n",
              "      <td>2.285</td>\n",
              "      <td>0</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>Union</th>\n",
              "      <td>26.930256</td>\n",
              "      <td>45087</td>\n",
              "      <td>2.931</td>\n",
              "      <td>2</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>Hampton</th>\n",
              "      <td>26.939454</td>\n",
              "      <td>45049</td>\n",
              "      <td>2.818</td>\n",
              "      <td>0</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>Jasper</th>\n",
              "      <td>28.094817</td>\n",
              "      <td>45053</td>\n",
              "      <td>3.006</td>\n",
              "      <td>2</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>Bamberg</th>\n",
              "      <td>30.703649</td>\n",
              "      <td>45009</td>\n",
              "      <td>2.902</td>\n",
              "      <td>0</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>Dillon</th>\n",
              "      <td>32.766392</td>\n",
              "      <td>45033</td>\n",
              "      <td>2.662</td>\n",
              "      <td>0</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>Kershaw</th>\n",
              "      <td>34.424195</td>\n",
              "      <td>45055</td>\n",
              "      <td>3.278</td>\n",
              "      <td>2</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>Calhoun</th>\n",
              "      <td>36.694939</td>\n",
              "      <td>45017</td>\n",
              "      <td>3.192</td>\n",
              "      <td>2</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>Orangeburg</th>\n",
              "      <td>38.795537</td>\n",
              "      <td>45075</td>\n",
              "      <td>2.901</td>\n",
              "      <td>0</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>Laurens</th>\n",
              "      <td>39.624397</td>\n",
              "      <td>45059</td>\n",
              "      <td>2.792</td>\n",
              "      <td>0</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>Allendale</th>\n",
              "      <td>40.158076</td>\n",
              "      <td>45005</td>\n",
              "      <td>2.431</td>\n",
              "      <td>0</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>Darlington</th>\n",
              "      <td>42.719819</td>\n",
              "      <td>45031</td>\n",
              "      <td>2.812</td>\n",
              "      <td>0</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>Georgetown</th>\n",
              "      <td>45.260607</td>\n",
              "      <td>45043</td>\n",
              "      <td>3.445</td>\n",
              "      <td>2</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>Oconee</th>\n",
              "      <td>48.786390</td>\n",
              "      <td>45073</td>\n",
              "      <td>2.818</td>\n",
              "      <td>0</td>\n",
              "    </tr>\n",
              "  </tbody>\n",
              "</table>\n",
              "</div>"
            ],
            "text/plain": [
              "              Euclidean Distance from Williamsburg COUNTYFP  Usage Index  Group\n",
              "Williamsburg                              0.000000    45089        3.295      2\n",
              "Colleton                                 16.709852    45029        3.129      2\n",
              "Saluda                                   17.057280    45081        3.063      2\n",
              "Clarendon                                17.925122    45027        2.818      0\n",
              "Lee                                      18.898024    45061        2.863      0\n",
              "Fairfield                                18.935085    45039        3.318      2\n",
              "Chesterfield                             20.095427    45025        2.866      0\n",
              "Marlboro                                 21.920866    45069        2.840      0\n",
              "Chester                                  23.574235    45023        2.865      0\n",
              "Edgefield                                24.201982    45037        3.021      2\n",
              "Barnwell                                 24.931590    45011        2.768      0\n",
              "Abbeville                                25.399076    45001        2.570      0\n",
              "Newberry                                 25.962924    45071        2.974      2\n",
              "Marion                                   26.588717    45067        2.285      0\n",
              "Union                                    26.930256    45087        2.931      2\n",
              "Hampton                                  26.939454    45049        2.818      0\n",
              "Jasper                                   28.094817    45053        3.006      2\n",
              "Bamberg                                  30.703649    45009        2.902      0\n",
              "Dillon                                   32.766392    45033        2.662      0\n",
              "Kershaw                                  34.424195    45055        3.278      2\n",
              "Calhoun                                  36.694939    45017        3.192      2\n",
              "Orangeburg                               38.795537    45075        2.901      0\n",
              "Laurens                                  39.624397    45059        2.792      0\n",
              "Allendale                                40.158076    45005        2.431      0\n",
              "Darlington                               42.719819    45031        2.812      0\n",
              "Georgetown                               45.260607    45043        3.445      2\n",
              "Oconee                                   48.786390    45073        2.818      0"
            ]
          },
          "metadata": {
            "tags": []
          },
          "execution_count": 413
        }
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "gVuYCnpFFx1z",
        "colab_type": "code",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 1000
        },
        "outputId": "01940b5e-b3e6-4796-f43b-f8d38c682322"
      },
      "source": [
        "closeness_scores.loc[donor_pool]\n",
        "#extremes\n",
        "X.loc[donor_pool+[county_of_interest]]\n",
        "state_counties_mask_usage.reset_index().set_index(\"Place\").loc[donor_pool]\n",
        "X.loc[[\"Aiken\", \"Sumter\", \"Beaufort\"]]\n",
        "#donor_pool = [\"Aiken\", \"Sumter\"]\n",
        "#mean_score\n",
        "\n",
        "#state_counties_mask_usage.sort_values(by = \"Usage Index\")\n",
        "donor_pool\n",
        "X.loc[donor_pool+[county_of_interest]]\n",
        "\n",
        "#donor_pool = [\"Oconee\"]\n",
        "#X[abs(X[\"Scaled magnitude of cases (1 day prior)\"]-4.045211)<1]\n",
        "state_counties_mask_usage.sort_values(by = \"Usage Index\")\n",
        "X.loc[donor_pool+[county_of_interest]]\n",
        "\n",
        "\n",
        "#closeness_scores2.join(state_counties_mask_usage.set_index(\"Place\")).sort_values(by = \"Usage Index\").head(15)\n",
        "\n",
        "state_counties_mask_usage.reset_index().set_index(\"Place\")"
      ],
      "execution_count": 443,
      "outputs": [
        {
          "output_type": "execute_result",
          "data": {
            "text/html": [
              "<div>\n",
              "<style scoped>\n",
              "    .dataframe tbody tr th:only-of-type {\n",
              "        vertical-align: middle;\n",
              "    }\n",
              "\n",
              "    .dataframe tbody tr th {\n",
              "        vertical-align: top;\n",
              "    }\n",
              "\n",
              "    .dataframe thead th {\n",
              "        text-align: right;\n",
              "    }\n",
              "</style>\n",
              "<table border=\"1\" class=\"dataframe\">\n",
              "  <thead>\n",
              "    <tr style=\"text-align: right;\">\n",
              "      <th></th>\n",
              "      <th>COUNTYFP</th>\n",
              "      <th>Usage Index</th>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>Place</th>\n",
              "      <th></th>\n",
              "      <th></th>\n",
              "    </tr>\n",
              "  </thead>\n",
              "  <tbody>\n",
              "    <tr>\n",
              "      <th>Abbeville</th>\n",
              "      <td>45001</td>\n",
              "      <td>2.570</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>Aiken</th>\n",
              "      <td>45003</td>\n",
              "      <td>2.917</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>Allendale</th>\n",
              "      <td>45005</td>\n",
              "      <td>2.431</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>Anderson</th>\n",
              "      <td>45007</td>\n",
              "      <td>2.993</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>Bamberg</th>\n",
              "      <td>45009</td>\n",
              "      <td>2.902</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>Barnwell</th>\n",
              "      <td>45011</td>\n",
              "      <td>2.768</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>Beaufort</th>\n",
              "      <td>45013</td>\n",
              "      <td>3.510</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>Berkeley</th>\n",
              "      <td>45015</td>\n",
              "      <td>3.189</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>Calhoun</th>\n",
              "      <td>45017</td>\n",
              "      <td>3.192</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>Charleston</th>\n",
              "      <td>45019</td>\n",
              "      <td>3.369</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>Cherokee</th>\n",
              "      <td>45021</td>\n",
              "      <td>2.974</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>Chester</th>\n",
              "      <td>45023</td>\n",
              "      <td>2.865</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>Chesterfield</th>\n",
              "      <td>45025</td>\n",
              "      <td>2.866</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>Clarendon</th>\n",
              "      <td>45027</td>\n",
              "      <td>2.818</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>Colleton</th>\n",
              "      <td>45029</td>\n",
              "      <td>3.129</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>Darlington</th>\n",
              "      <td>45031</td>\n",
              "      <td>2.812</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>Dillon</th>\n",
              "      <td>45033</td>\n",
              "      <td>2.662</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>Dorchester</th>\n",
              "      <td>45035</td>\n",
              "      <td>3.101</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>Edgefield</th>\n",
              "      <td>45037</td>\n",
              "      <td>3.021</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>Fairfield</th>\n",
              "      <td>45039</td>\n",
              "      <td>3.318</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>Florence</th>\n",
              "      <td>45041</td>\n",
              "      <td>2.929</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>Georgetown</th>\n",
              "      <td>45043</td>\n",
              "      <td>3.445</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>Greenville</th>\n",
              "      <td>45045</td>\n",
              "      <td>3.055</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>Greenwood</th>\n",
              "      <td>45047</td>\n",
              "      <td>2.844</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>Hampton</th>\n",
              "      <td>45049</td>\n",
              "      <td>2.818</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>Horry</th>\n",
              "      <td>45051</td>\n",
              "      <td>3.148</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>Jasper</th>\n",
              "      <td>45053</td>\n",
              "      <td>3.006</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>Kershaw</th>\n",
              "      <td>45055</td>\n",
              "      <td>3.278</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>Lancaster</th>\n",
              "      <td>45057</td>\n",
              "      <td>3.212</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>Laurens</th>\n",
              "      <td>45059</td>\n",
              "      <td>2.792</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>Lee</th>\n",
              "      <td>45061</td>\n",
              "      <td>2.863</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>Lexington</th>\n",
              "      <td>45063</td>\n",
              "      <td>3.242</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>McCormick</th>\n",
              "      <td>45065</td>\n",
              "      <td>2.953</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>Marion</th>\n",
              "      <td>45067</td>\n",
              "      <td>2.285</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>Marlboro</th>\n",
              "      <td>45069</td>\n",
              "      <td>2.840</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>Newberry</th>\n",
              "      <td>45071</td>\n",
              "      <td>2.974</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>Oconee</th>\n",
              "      <td>45073</td>\n",
              "      <td>2.818</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>Orangeburg</th>\n",
              "      <td>45075</td>\n",
              "      <td>2.901</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>Pickens</th>\n",
              "      <td>45077</td>\n",
              "      <td>3.031</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>Richland</th>\n",
              "      <td>45079</td>\n",
              "      <td>3.399</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>Saluda</th>\n",
              "      <td>45081</td>\n",
              "      <td>3.063</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>Spartanburg</th>\n",
              "      <td>45083</td>\n",
              "      <td>3.008</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>Sumter</th>\n",
              "      <td>45085</td>\n",
              "      <td>2.840</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>Union</th>\n",
              "      <td>45087</td>\n",
              "      <td>2.931</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>Williamsburg</th>\n",
              "      <td>45089</td>\n",
              "      <td>3.295</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>York</th>\n",
              "      <td>45091</td>\n",
              "      <td>3.218</td>\n",
              "    </tr>\n",
              "  </tbody>\n",
              "</table>\n",
              "</div>"
            ],
            "text/plain": [
              "             COUNTYFP  Usage Index\n",
              "Place                             \n",
              "Abbeville       45001        2.570\n",
              "Aiken           45003        2.917\n",
              "Allendale       45005        2.431\n",
              "Anderson        45007        2.993\n",
              "Bamberg         45009        2.902\n",
              "Barnwell        45011        2.768\n",
              "Beaufort        45013        3.510\n",
              "Berkeley        45015        3.189\n",
              "Calhoun         45017        3.192\n",
              "Charleston      45019        3.369\n",
              "Cherokee        45021        2.974\n",
              "Chester         45023        2.865\n",
              "Chesterfield    45025        2.866\n",
              "Clarendon       45027        2.818\n",
              "Colleton        45029        3.129\n",
              "Darlington      45031        2.812\n",
              "Dillon          45033        2.662\n",
              "Dorchester      45035        3.101\n",
              "Edgefield       45037        3.021\n",
              "Fairfield       45039        3.318\n",
              "Florence        45041        2.929\n",
              "Georgetown      45043        3.445\n",
              "Greenville      45045        3.055\n",
              "Greenwood       45047        2.844\n",
              "Hampton         45049        2.818\n",
              "Horry           45051        3.148\n",
              "Jasper          45053        3.006\n",
              "Kershaw         45055        3.278\n",
              "Lancaster       45057        3.212\n",
              "Laurens         45059        2.792\n",
              "Lee             45061        2.863\n",
              "Lexington       45063        3.242\n",
              "McCormick       45065        2.953\n",
              "Marion          45067        2.285\n",
              "Marlboro        45069        2.840\n",
              "Newberry        45071        2.974\n",
              "Oconee          45073        2.818\n",
              "Orangeburg      45075        2.901\n",
              "Pickens         45077        3.031\n",
              "Richland        45079        3.399\n",
              "Saluda          45081        3.063\n",
              "Spartanburg     45083        3.008\n",
              "Sumter          45085        2.840\n",
              "Union           45087        2.931\n",
              "Williamsburg    45089        3.295\n",
              "York            45091        3.218"
            ]
          },
          "metadata": {
            "tags": []
          },
          "execution_count": 443
        }
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "e6DuDJHfSqph",
        "colab_type": "code",
        "colab": {}
      },
      "source": [
        "closeness_scores2.join(state_counties_mask_usage.set_index(\"Place\")).sort_values(by = closeness_scores2.columns[0]).head(15)\n",
        "\n",
        "donor_pool = [\"Lee\", \"Marlboro\", \"Clarendon\"]"
      ],
      "execution_count": 438,
      "outputs": []
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "qvCQ6zlpeXWz",
        "colab_type": "code",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 204
        },
        "outputId": "4cfbc94e-8602-4350-a6cb-b8e7bf26f39f"
      },
      "source": [
        "state = pd.read_csv(URL1, index_col=0)\n",
        "looking_for_alternative = \"cases\"\n",
        "state=state[state['state']==target].pivot_table(index='county', columns='date', values= looking_for_alternative).fillna(0)\n",
        "data = state[[i for i in state.columns if \n",
        "              (dtm.strptime(i, format)<=dtm.strptime(end, format)) and \n",
        "              (dtm.strptime(i, format)>=dtm.strptime(start, format))]].loc[donor_pool + [county_of_interest]]\n",
        "data[[\"2020-05-30\", \"2020-05-31\", \"2020-06-01\", \"2020-06-02\"]]"
      ],
      "execution_count": 439,
      "outputs": [
        {
          "output_type": "execute_result",
          "data": {
            "text/html": [
              "<div>\n",
              "<style scoped>\n",
              "    .dataframe tbody tr th:only-of-type {\n",
              "        vertical-align: middle;\n",
              "    }\n",
              "\n",
              "    .dataframe tbody tr th {\n",
              "        vertical-align: top;\n",
              "    }\n",
              "\n",
              "    .dataframe thead th {\n",
              "        text-align: right;\n",
              "    }\n",
              "</style>\n",
              "<table border=\"1\" class=\"dataframe\">\n",
              "  <thead>\n",
              "    <tr style=\"text-align: right;\">\n",
              "      <th>date</th>\n",
              "      <th>2020-05-30</th>\n",
              "      <th>2020-05-31</th>\n",
              "      <th>2020-06-01</th>\n",
              "      <th>2020-06-02</th>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>county</th>\n",
              "      <th></th>\n",
              "      <th></th>\n",
              "      <th></th>\n",
              "      <th></th>\n",
              "    </tr>\n",
              "  </thead>\n",
              "  <tbody>\n",
              "    <tr>\n",
              "      <th>Lee</th>\n",
              "      <td>217.0</td>\n",
              "      <td>218.0</td>\n",
              "      <td>218.0</td>\n",
              "      <td>221.0</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>Marlboro</th>\n",
              "      <td>161.0</td>\n",
              "      <td>171.0</td>\n",
              "      <td>173.0</td>\n",
              "      <td>183.0</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>Clarendon</th>\n",
              "      <td>301.0</td>\n",
              "      <td>301.0</td>\n",
              "      <td>301.0</td>\n",
              "      <td>302.0</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>Williamsburg</th>\n",
              "      <td>233.0</td>\n",
              "      <td>236.0</td>\n",
              "      <td>238.0</td>\n",
              "      <td>245.0</td>\n",
              "    </tr>\n",
              "  </tbody>\n",
              "</table>\n",
              "</div>"
            ],
            "text/plain": [
              "date          2020-05-30  2020-05-31  2020-06-01  2020-06-02\n",
              "county                                                      \n",
              "Lee                217.0       218.0       218.0       221.0\n",
              "Marlboro           161.0       171.0       173.0       183.0\n",
              "Clarendon          301.0       301.0       301.0       302.0\n",
              "Williamsburg       233.0       236.0       238.0       245.0"
            ]
          },
          "metadata": {
            "tags": []
          },
          "execution_count": 439
        }
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "_iTuZpfIeXW1",
        "colab_type": "text"
      },
      "source": [
        "## Performing Synthetic Control"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "eB8NzVJGeXW3",
        "colab_type": "code",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 419
        },
        "outputId": "3978ed6e-d724-46ae-c559-7ff75f409594"
      },
      "source": [
        "allColumns = data.columns.values\n",
        "years = list(data.columns)\n",
        "caStateKey = county\n",
        "otherStates = [i for i in data.index if i!=county_of_interest]\n",
        "p = 1.0\n",
        "yearTrainEnd = intervention\n",
        "#*******************************************************************************\n",
        "testYears=list(years)[list(years).index(date):len(list(years))]\n",
        "trainingYears=list(years)[list(years).index(start):list(years).index(date)]\n",
        "trainDF = data[trainingYears].T.reset_index(drop = True)\n",
        "testDF = data[testYears].T.reset_index(drop = True)\n",
        "trainDF = data[trainingYears].T.reset_index(drop = True)\n",
        "testDF = data[testYears].T.reset_index(drop = True)\n",
        "\n",
        "trainDF"
      ],
      "execution_count": 440,
      "outputs": [
        {
          "output_type": "execute_result",
          "data": {
            "text/html": [
              "<div>\n",
              "<style scoped>\n",
              "    .dataframe tbody tr th:only-of-type {\n",
              "        vertical-align: middle;\n",
              "    }\n",
              "\n",
              "    .dataframe tbody tr th {\n",
              "        vertical-align: top;\n",
              "    }\n",
              "\n",
              "    .dataframe thead th {\n",
              "        text-align: right;\n",
              "    }\n",
              "</style>\n",
              "<table border=\"1\" class=\"dataframe\">\n",
              "  <thead>\n",
              "    <tr style=\"text-align: right;\">\n",
              "      <th>county</th>\n",
              "      <th>Lee</th>\n",
              "      <th>Marlboro</th>\n",
              "      <th>Clarendon</th>\n",
              "      <th>Williamsburg</th>\n",
              "    </tr>\n",
              "  </thead>\n",
              "  <tbody>\n",
              "    <tr>\n",
              "      <th>0</th>\n",
              "      <td>0.0</td>\n",
              "      <td>0.0</td>\n",
              "      <td>0.0</td>\n",
              "      <td>0.0</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>1</th>\n",
              "      <td>0.0</td>\n",
              "      <td>0.0</td>\n",
              "      <td>0.0</td>\n",
              "      <td>0.0</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>2</th>\n",
              "      <td>0.0</td>\n",
              "      <td>0.0</td>\n",
              "      <td>0.0</td>\n",
              "      <td>0.0</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>3</th>\n",
              "      <td>0.0</td>\n",
              "      <td>0.0</td>\n",
              "      <td>0.0</td>\n",
              "      <td>0.0</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>4</th>\n",
              "      <td>0.0</td>\n",
              "      <td>0.0</td>\n",
              "      <td>0.0</td>\n",
              "      <td>0.0</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>...</th>\n",
              "      <td>...</td>\n",
              "      <td>...</td>\n",
              "      <td>...</td>\n",
              "      <td>...</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>79</th>\n",
              "      <td>209.0</td>\n",
              "      <td>137.0</td>\n",
              "      <td>296.0</td>\n",
              "      <td>227.0</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>80</th>\n",
              "      <td>212.0</td>\n",
              "      <td>149.0</td>\n",
              "      <td>296.0</td>\n",
              "      <td>229.0</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>81</th>\n",
              "      <td>217.0</td>\n",
              "      <td>161.0</td>\n",
              "      <td>301.0</td>\n",
              "      <td>233.0</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>82</th>\n",
              "      <td>218.0</td>\n",
              "      <td>171.0</td>\n",
              "      <td>301.0</td>\n",
              "      <td>236.0</td>\n",
              "    </tr>\n",
              "    <tr>\n",
              "      <th>83</th>\n",
              "      <td>218.0</td>\n",
              "      <td>173.0</td>\n",
              "      <td>301.0</td>\n",
              "      <td>238.0</td>\n",
              "    </tr>\n",
              "  </tbody>\n",
              "</table>\n",
              "<p>84 rows × 4 columns</p>\n",
              "</div>"
            ],
            "text/plain": [
              "county    Lee  Marlboro  Clarendon  Williamsburg\n",
              "0         0.0       0.0        0.0           0.0\n",
              "1         0.0       0.0        0.0           0.0\n",
              "2         0.0       0.0        0.0           0.0\n",
              "3         0.0       0.0        0.0           0.0\n",
              "4         0.0       0.0        0.0           0.0\n",
              "..        ...       ...        ...           ...\n",
              "79      209.0     137.0      296.0         227.0\n",
              "80      212.0     149.0      296.0         229.0\n",
              "81      217.0     161.0      301.0         233.0\n",
              "82      218.0     171.0      301.0         236.0\n",
              "83      218.0     173.0      301.0         238.0\n",
              "\n",
              "[84 rows x 4 columns]"
            ]
          },
          "metadata": {
            "tags": []
          },
          "execution_count": 440
        }
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "uf6x0mWLeXW9",
        "colab_type": "code",
        "colab": {}
      },
      "source": [
        "singvals = 4\n",
        "rscModel = RobustSyntheticControl(caStateKey, singvals, len(trainDF), probObservation=1.0, modelType='svd', svdMethod='numpy', otherSeriesKeysArray=otherStates)\n",
        "rscModel.fit(trainDF)\n",
        "denoisedDF = rscModel.model.denoisedDF()\n",
        "predictions = []\n",
        "predictions = np.dot(testDF[otherStates], rscModel.model.weights)\n",
        "actual = data.loc[caStateKey]\n",
        "model_fit = np.dot(trainDF[otherStates][:], rscModel.model.weights)"
      ],
      "execution_count": 441,
      "outputs": []
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "LgYLhVcDeXXA",
        "colab_type": "text"
      },
      "source": [
        "## Plotting graph"
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "P4SU9VeDeXXA",
        "colab_type": "code",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 638
        },
        "outputId": "0395210f-0851-4b0d-bb21-9faacb24e200"
      },
      "source": [
        "# fig,ax = plt.subplots(1,1)\n",
        "# fig.autofmt_xdate()\n",
        "plt.figure(figsize=(15,10))\n",
        "tick_spacing = 1\n",
        "intDate=yearTrainEnd\n",
        "# this is a bug in matplotlib\n",
        "label_markings = np.insert(years[::tick_spacing], 0, 'dummy')\n",
        "\n",
        "# ax.set_xticks(np.arange(len(label_markings)))\n",
        "#ax.set_xticklabels(label_markings, rotation=45)\n",
        "# ax.xaxis.set_major_locator(ticker.MultipleLocator(tick_spacing))\n",
        "\n",
        "plt.plot(years, actual ,label='Actual County: '+county_of_interest)\n",
        "\n",
        "#for i in donor_pool:\n",
        "#  plt.plot(years, data.loc[i], label = i)\n",
        "plt.xlabel('Time Period')\n",
        "plt.ylabel('No. of confirmed Cases')\n",
        "plt.plot(trainingYears, model_fit, label='fitted model')\n",
        "plt.plot(testYears, predictions, label='counterfactual')\n",
        "plt.title(caStateKey+', Singular Values used: '+str(singvals))\n",
        "\n",
        "xposition = pd.to_datetime(yearTrainEnd,  errors='coerce')\n",
        "plt.axvline(x=str(intDate), color='k', linestyle='--', linewidth=4)\n",
        "plt.grid()\n",
        "plt.legend()"
      ],
      "execution_count": 442,
      "outputs": [
        {
          "output_type": "execute_result",
          "data": {
            "text/plain": [
              "<matplotlib.legend.Legend at 0x7f84c93b0668>"
            ]
          },
          "metadata": {
            "tags": []
          },
          "execution_count": 442
        },
        {
          "output_type": "display_data",
          "data": {
            "image/png": "iVBORw0KGgoAAAANSUhEUgAAA3sAAAJcCAYAAABAE73ZAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4yLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+WH4yJAAAgAElEQVR4nOzdd3hUZfrG8e+bQkIJhBoElNBBSQiEZigGEUVBbKACi2ABdFdwF11dy+riomtFECsoqCyCiisqiGIh9BpRpNfQIiVAIBWSzPv7Y4b5JUgJJTmTyf25rlzOnPbcZyZe5vF9zznGWouIiIiIiIj4lwCnA4iIiIiIiMjFp2ZPRERERETED6nZExERERER8UNq9kRERERERPyQmj0RERERERE/pGZPRERERETED6nZExEpIYwxnYwxG/O9TzLGXON5/S9jzH89ry8zxqQbYwKLMVukMcYaY4KKq+b5Msb0N8bMKYY6PvOZ5P/98CfGmHhjzG6nc4iI+Co1eyIiDjHGPG6MmX3Sss2nWXantXaBtbbJ2Y5rrd1pra1grc272JlLCmNMR2PMYmPMEWPMIWPMImNMGwBr7RRr7bVOZzwXxpj2xpgMY0yFU6xbZYx50IlcJZUx5ipPIz7K6SwiIkVJzZ6IiHPmA3EnRuCMMZcAwUDLk5Y19GxbKlzoSJgxpiIwExgHVAFqAyOBYxeernic/BlYa5cCu4HeJ23XHLgcmFp86Uo2Y0wwMBZY5nQWEZGipmZPRMQ5K3A3dzGe952AucDGk5ZttdYmF3bK2snTB40xdxtj1htj0owx24wxQ/NtG2+M2W2MedQYs98Y87sx5mZjzA3GmE2eUbEn8m3f1hiz0hhz1Bizzxgz+qTy9xhjkj3HeSTffh/kH0U5+Vw8U1IfM8asBjKMMUHGmLuMMTuMMQeNMf/MP231LBoDWGunWmvzrLVZ1to51trVnlqDjDEL89W2xpj7PSOoqcaYN40xxrMu0BjzqjEmxRiz3Rjz4EmfbYFMZ5ouWcjv4TFjzF5g0ikO8SFw10nL7gK+sdYeNMaMNcbs8nw3icaYTqfJ8YffI1NwSnCAMeYfxpitns/+U2NMFc+6UGPMfz3LU40xK4wxEaepY40xDfO99/4OGGOqGWNmeo5xyBizwBgT4FlXyxjzuTHmgOczH57vGGU9xzlsjFkHtDlV7bN4GJgDbDiPfUVEShQ1eyIiDrHWHsc9utDZs6gzsABYeNKyCx3V2w/0BCoCdwOvGWNa5VtfEwjFPQL2NDAB+BMQi7vZ/Kcxpp5n27HAWGttRaAB8OlJtboAjYBrgccK2Zyd0BfoAYTjbtjeAvoDlwCVPPkKYxOQZ4z50BhzvTGmciH26Ym7cYgGbgeu8ywfDFyPu/luBdxcyAynUpjvoQpQFxhyiv0nA52NMZeCuykD+uFuAsH9Pw9iPMf4GPjMGBN6HjmH4T7Pq4BawGHgTc+6gbi/i0uBqsD9QNZ51HgY90hldSACeAKwnnP6GvgV9/fdFfirMebE9/EM7t+7Bri/o4H5D2qMecsY89bpihpj6gL3AM+eR2YRkRJHzZ6IiLPm8f+NXSfczd6Ck5bNu5AC1tpZ1tqt1m0e7lGN/KM+OcBz1tocYBpQDXdDl2atXQusA1rk27ahMaaatTbdM70wv5HW2gxr7W+4R6f6nkPU1621u6y1WbinK35trV3oaYqfBmwhz/co0NGz/QTggDHmq9ONQHm8YK1NtdbuxD26emJk9Xbcn8Vua+1h4IVzOJ+Tc53te3ABz1hrj3k+g5P33wUkAAM8i7oCIcAsz/r/WmsPWmtzrbWvetad9RrPU7gfeNJzzseAfwG9PaOZObibvIaeUdNEz+d9rnJwN/F1rbU5nutRLe6Gu7q19llr7XFr7Tbc3+Gdnv1ux/27esjzebye/6DW2j9ba/98hrqvA/+01qafR2YRkRJHzZ6IiLPmAx090+SqW2s3A4txX8tXBWjOBY7seUa3lnqmy6UCN+Bu6E44mO9mLieajH351mcBJ24Mci/uUbcNnil8PU8qtyvf6x24R4YKK/++tfK/t9ZmAgcLeyBr7Xpr7SBrbR3cn2EtYMwZdtmb73Um/3++BXKc9PqcFOJ7OGCtzT7LYT7k/5u9AcA0T5OOMeYRzzTRI57jVzrp+IVVF/jCM8UyFVgP5OEegZsMfAdM80zXfcm4r4E7Vy8DW4A5nimt/8hXu9aJ2p76T3hqwx+/jx2FLWiMuREIs9Z+ch55RURKJDV7IiLOWoL7j/LBwCLwjkwle5YlW2u3n+/BjTEhwOfAK0CEtTYc+AYw53M8a+1ma21foAbwIjDdGFM+3yaX5nt9Ge7zAMgAyuVbV/NUh8/3+negTr7zKIt7ROl8Mm8APsDd9J2rAjkoeH5QuPMq7PdQmJHL/wF1jDFdgFvxTOH0XJ/3KO6Rr8qe4x/h1N9zgczGfTOg6vnW7wKut9aG5/sJtdbu8YzCjbTWXg7E4Z6WevJ1hCdkcprPxjNq/LC1tj7QCxhhjOnqqb39pNph1tobPLv+zh9/xwqrK9DaGLPXc13kHbiniH55DscQESlR1OyJiDjIM11vJTAC9/TNExZ6ll3o9XplcE/nOwDkGmOux3093XkxxvzJGFPdWusCUj2LXfk2+acxppwx5grc16WdGEX5BbjBGFPFGFMT+OtZSk0HbjTGxBljyuCeSuhtXDw3GTllc2SMaWqMedgYU8fz/lLc00lPnnJaGJ8CDxljahtjwoHHTlr/C3CnMSbYGNOak+6Wmc9F+R6stRm4P5tJwA5r7UrPqjAg13P8IGPM07ivDTyVTUCoMaaHZ1TuKU+2E94BnvNc34Yxprox5ibP6y7GmChPg3gU93RMF6f2C9DPuG9y0x33NYB4jtPTGNPQGGNwN6V5nuMsB9KM+0Y1ZT37Njeex2bg/j4eN8ZU9ny/wwrzuXn8E/eodIzn5yvcU0TvPodjiIiUKGr2REScNw/3SNnCfMsWeJZdULNnrU0DhuP+I/kw7ht6fHUBh+wOrDXGpOO+WcudJ11fNg/39LwfgVestSceXj4Z9003knBfq3bGqXSeawWH4b6G8HcgHfcNTk48PuFS3NNdTyUNaAcsM8Zk4G7y1uC+Kci5muDJuxpYhXs0Lhd3cwLuBqIB7s92JO4bo5zqfC7m9/Ah7umOH+Vb9h3wLe5GbgeQzWmmnFprjwB/Bt4D9uAe6ct/d86xnmxzjDFpuD+/dp51NXE3m0dxT++ch/u7PZWHgBtx/0+B/sCMfOsaAT/g/l6XAG9Za+d6phP3xN2MbQdSPDkrefYb6Tm/7bi/lwK1jTHvGGPeOc15p1lr9574wT09OcNae+g0+UVESjzjvh5aRETEdxn3w8RTgUbW2u3GmPeAz6y13xVzjuuBd6y1dYuzroiIyPlQsyciIj7Jc0ONH3FP33wV9+hSK1uM/+HyXCvYBfcoUgTu6+6WWmvPNg1VRETEcZrGKSIivuom3Dd4ScY97e/O4mz0PAzuqYOHcU/jXI/7MRAiIiI+TyN7IiIiIiIifkgjeyIiIiIiIn4oyOkAF6JatWo2MjLS6Rh/kJGRQfny5Qu9/Hz28cUa/n5+qlEyaqtGyaitGr5Vw9/Pz9dqJCYmFlgfGxtbbLVVwz9q+Pv5lbQaTktMTEyx1lY/5UprbYn9iY2Ntb5o7ty557T8fPbxxRr+fn6qUTJqq0bJqK0avlXD38/P12oABX6Ks7Zq+EcNfz+/klbDacBKe5p+SdM4RURERERE/JCaPRERERERET+kZk9ERERERMQPlegbtJxKTk4Ou3fvJjs727EMlSpVYv369YVefj77+GINfz+/4qgRGhqKMeaUtUVEREREzoXfNXu7d+8mLCyMyMhIx/5oTktLIywsrNDLz2cfX6zh7+dX1DWstRw8eNAn7/IkIiIiIiWP303jzM7OpmrVqhodkRLHGEPVqlUJDAx0OoqIiIiI+AG/a/YANXpSYul3V0REREQuFr9s9kREREREREo7NXtFZMaMGRhj2LBhw1m3HTNmDJmZmedd64MPPuDBBx885brZs2fTunVrLr/8cjp27MjDDz983nVO53zzW2uJjIzk8OHDAPz+++8YY1i4cKF3m+rVq3Pw4EHuu+8+1q1bB0BkZCQpKSkAVKhQAYDk5GQGDBhwoadyWmf6jEVEREREfJGavSIydepUOnbsyNSpU8+67ZgxY8jKyrroGdasWcODDz7If//7X9atW8e8efNo2LDhRa9zvs2eMYY2bdqwZMkSABYvXkzLli1ZvHgxABs3bqRq1apUrVqV9957j8svv/y0x6pVqxaTJ08+vxMoBrm5uU5HEBEREZFSRs1eEUhPT2fhwoW8//77TJs2zbs8Ly+PRx55hObNmxMdHc24ceN4/fXXSU5OpkePHnTp0gX4/9EqcI8QDho0CICvv/6adu3a0bJlS3r16sW+ffvOmOOll17iySefpGnTpgAEBgbywAMPAJCUlMTVV19NdHQ0N954Izt37gRg0KBBTJ8+3XuME1kSEhKIj4+nd+/exMbG0r9/f6y13vxdunShR48eTJw4kb/+9a/e/SdMmMDf/va302Zs166dt7lbvHgxf/vb3wo0fx06dAAgPj6elStXnvY4SUlJtGvXzvu6U6dOtGrVilatWrFs2TLvOVx11VXcdNNN1K9fn2eeeYYpU6bQtm1boqKi2Lp1KwBffPEFzZs3p0WLFnTu3NlbY9euXdxwww00atSIkSNHems1b97cu80rr7zCv/71L2/mxx57jNatWzN27FhWrFhBdHQ0MTExPPXUUwX2ExERERG52Pzu0Qv5jfx6LeuSj17UY15eqyLP3HjFGbeZNWsW3bt3p3HjxlStWpXExERiY2OZNGkSSUlJ/PLLLwQFBXHo0CGqVKnC6NGjmTVrFpGRkWc8bseOHVm6dCnGGN544w1eeuklXn311dNuv2bNmtNO2xw2bBgDBw5k4MCBvPXWWwwfPpwZM2acsf6qVatYu3YtYWFhdO/enUWLFjF8+HBGjx7N3LlzCQkJwRjDc889x8svv0xwcDCTJk3i3XffZdSoUXTo0IFevXoVOGb79u15+eWXAVi+fDkjR45k7NixgLvZi4uLO2OmU6lRowbff/89oaGhbN68mTvuuIOff/4ZgF9//ZX169dTpUoV6tWrR1hYGMuXL2fs2LGMGzeOMWPG8OKLLzJnzhxq165Namqq97jLly9nyZIlRERE0KZNG3r06EG1atXOmOX48ePeJrV58+ZMmDCBK6+8khEjRpzzeYmIiIiInAuN7BWB6dOnc+eddwJw5513eqdyJiQkMHToUIKC3D12lSpVzum4u3fv5rrrriMqKoqxY8eydu3a8864ZMkS+vXr582Y/zq502nbti116tQhICCAmJgYkpKS/rBNhQoVuPrqq5k5cyYbNmwgJyeHqKgonnrqqT80egCtWrVi1apVZGRkkJOTQ4UKFahfvz5bt24tMLJ3LnJychg8eDBRUVH06dOnwHWTbdq04ZJLLiEkJIR69epx7bXXAhAVFeU9n/bt2zNo0CAmTJhAXl6ed99u3bpRtWpVypYty6233lqoz+y2224DIDU1lbS0NK688koA+vTpc87nJSIi/sFai7WWuXPnYq11Oo6I+DG/Htk72whcUTh06BDz589n/fr1GGPIy8vDGOMdvSqM/Lffz87O9r4eNmwYI0aMoFevXnzzzTe89NJLZzzOFVdcQWJiIi1atCh07aCgIFwuFwAul4vjx49714WEhHhfBwYGnvY6tPvuu4/nn3+epk2bcvfdd5+xXrly5WjUqBETJ06kVatWgLvZmjNnDvv376dJkyaFzn7Ca6+9RkREBL/++isul4vQ0NBTnkNAQID3fUBAgPd8xowZw7p165g1axaxsbEkJiYCf3wsgjGmwOcFBb+vE+cnIiIiIuIEjexdZCdG9Xbs2EFSUhK7du2iXr16LFiwgC5duvDuu+96m4pDhw4BEBYWRlpamvcYERERrF+/HpfLxcyZM73Ljxw5Qu3atQH4+OOPz5rl73//O88//zybNm0C3M3bO++8A0BcXJz3esJPP/2UTp06Ae47XZ5obr755htycnLOWufk/O3atWPXrl18/PHH9O3b96z7x8XFMWbMGO+o15VXXsnbb79N+/btz+u5c0eOHOGSSy4hICCAyZMnFxidK4xt27bRrl07nn32WapXr86uXbsA+P777zl06BBZWVnMmDGDDh06EBERwf79+zl48CDHjh0r8H3lFx4eTlhYmPf6wc8///ycz0tERERE5Fyo2bvIpk6dSs+ePQssu+2225g6dSoDBw7ksssuIzo6mhYtWngbtiFDhnDrrbd6b9Dywgsv0LNnT+Li4oiIiPAe51//+hd9+vQhNjaWqlWrnjVLdHQ0Y8aMoW/fvjRr1ox27dqxbds2AMaNG8ekSZOIjo5m2rRp3uvkBg8ezLx582jRogXLly+nfPnyZ60zZMgQunfvTo8ePbzLbr/9djp06EDlypUBGDVqFF999dUp9+/QoQPbtm3zNnutWrUiOTn5vK7XA/jzn//Mhx9+SIsWLdiwYUOhziG/f/7zn0RFRdG8eXPi4uK8I6Nt27ZlwIABREdHc9ttt9G6dWuCg4N5+umn6dKlC926dfPeDOdU3n//fQYPHkxMTAwZGRlUqlTpvM5PRERERKQw/HoapxPmzp1bYJQLYPjw4QCkpaUxevRoRo8eXWD9sGHDGDRoEGFhYQD07t2b3r17e/c5sfymm27ipptu+sPyQYMGee/YebKePXt6m8/8+9StW5effvrpD8sjIiJYunSpd/lrr70GuO8sGR8f7z3uG2+8USD/sGHDCpz3woULC9yF86mnnvLWOFmfPn0KXLMQEhJCSkpKge0TEhK8mfJfK5ieng64RyRPjJo1atSI1atXF6h9qnP45ptvvDXyr5syZcofsp74jPN/VicMHz6cu++++w/LExISCnwmV1xxhTfXyJEjKVOmzCk/DxERERGRi0Eje3JRpaam0rhxY8qWLUvXrl2djuNTZs2aRUxMDM2bN2fx4sXeJlREREREpChoZE8uqvDwcO81glLQHXfcwR133AFwyhFCEREpHYYMGQJAcnIyH3/8MePHj3c4kYj4KzV7IiIiIsVowoQJBd6r2RORoqJpnCIiIiIiIn5IzZ6IiIiIiMhZ5Lrs2TfyMWr2RERERERETiM18zgPTVvF+NXHnI5yztTsFYG3336bZs2a0b9/f7766iteeOEFAGbOnMm6deu8233wwQckJyef07GTkpJo3rz5Rc17KhUqVLgo24iIiIiIlFQ/bdjHta/NZ9bq36ldIQBXCRvd0w1aisB7773HTz/9RJ06dQDo1asX4G72goODufzyywF3s9e8eXNq1arlWFYRERERESnoaHYO//56HZ8l7qZpzTAmDmpDyuZVBAQYp6OdEzV7F9n9999PUlIS119/Pffccw+VK1dm5cqV9OvXj2+++YbFixczatQo+vbty8qVK+nfvz9ly5Zlzpw5bNq0iREjRpCenk61atX44IMPqFChAomJidxzzz0AXHvttaesm5CQwDPPPEOFChVYv349t99+O1FRUYwdO5asrCxmzJhBjRo1SEpK4p577iElJYXq1aszadIkKleuzPbt2+nXrx/p6eneB7ef8PLLL/Ppp59y7NgxbrnlFh555JEi/xxFRERERJywJiWPJ16bz96j2fylSwOGd21ESFAgCZudTnbu/LvZm/0P2PvbxT1mzSi4/oXTrn7nnXeYPXs2c+fO9TZsAHFxcdxwww3ccsst9O7d2x1v9mxeeeUVWrduzaFDhxg2bBhffvkl1atX55NPPuHJJ59k7Nix3H333bzxxht07tyZv//976et/euvv7JixQouu+wy6tevz3333cfy5csZO3Ys48aN49///jfDhg1j4MCBDBw4kIkTJzJ8+HAmT57MQw89xAMPPMBdd93Fm2++6T3mnDlz2Lx5M8uXL8daS69evVi0aBHdu3e/OJ+niIiIiIgPyHNZ/j1zHR+szKZ+9fJ8/kAcLS+r7HSsC+LfzV4JsnnzZtasWUO3bt0AyMvL45JLLiE1NZXU1FQ6d+4MwIABA5g9e/Ypj9GmTRtq1qxJSEgIDRo08I4CRkVFMXfuXACWLFnC//73P++xHn30UQAWLVrE559/7l3+2GOPAe5mb86cObRs2RKA9PR0tm7dWhQfgYiIiIiII/Jclr9P/5X//byHbnWDGHdfJ0KDA52OdcH8u9k7wwicr7HWcsUVV7BkyZICy3ft2lXoY4SEhHhfBwQEeN8HBASQm5t71v2N+eMcZGstjz/+OEOHDvUuS0tLK3QmERERERFf5nJZ/vH5av738x5GdGtMdOAev2j0QHfjLFYVKlQo0CiFhYV53zdq1IgDBw54m72cnBzWrl1LeHg44eHhLFy4EIApU6ZcUIa4uDimTZvmPVanTp0A6NChQ4HlJ1x33XVMnDiR9PR0APbs2cOBAwcuKIOIiIiIiC9wuSxPfPEbnyXu5qGujRjetZHTkS4qNXvFqHfv3rz88su0bNmSrVu3MmjQIO6//35iYmLIy8tj+vTpPPbYY7Ro0YKYmBgWL14MwKRJk/jLX/5CTEwM1l7Y7V7HjRvHpEmTiI6OZvLkyYwdOxaAsWPH8uabbxIVFcWePXu821977bX069ePK6+8kqioKHr37q2RPREREREp8ay1PPXlGqat2MWDXRry12v8q9EDf5/G6ZA1a9YQFhYGwKBBgxg0aBAA7du3L/CcvQYNGnDbbbcB7qmRMTExzJ8/v8Cx0tLSiI2N5ddff/Uue+mll/7QcMXHxxMfH+9dnpCQcMp1devW5aeffvpDjXr16hWYQjpq1CjvsR566CEeeuihAtsD3tE+EREREZGSxFrL5PXH+WnnTh6Ib8DD1zY+5SVNJZ1G9kREREREpNSw1jLy63X8tDOXIZ3r8+h1Tfyy0QON7ImIiIgUq5UrV3r/2bp1a4fTiJQ+CZsO8MHiJLrVDeLx65v6baMHavZEREREilVsbCzw/5dqiEjx+mhxEtXDQrijSaBfN3qgaZwiIiIiIlJKJKVkkLDpAP3aXkZQgH83eqBmT0RERERESonJS3cQaAz92l3mdJRioWZPRERERET8XubxXD5duYvuzWsSUTHU6TjFQs1eCTRmzBgyMzPPeb9NmzYRExPjfc7fuUhISPA+9+98VahQ4YL2FxERERE5XzNWJZOWncuguEinoxQbNXsl0JgxY8jKyjqnffLy8pg5cya9e/dm1apVNGjQ4Jz2vxjNnoiIiMD48eMZP348X3/9NePHj3c6jkipYK3loyVJXH5JRWLrVnY6TrFRs1dEPvroI6Kjo2nRogUDBgwgKSmJnj17Eh0dTdeuXdm5cyfgfuj69OnTvfudGP1KSEggPj6eAQMG0LRpU/r374+1ltdff53k5GR69OhBly5dAJgzZw5XXnklrVq14q677vI+7DwyMpLHHnuMVq1a8cknn/DWW2/x9ttve/e7+eabiY2NpW3btgX+Y/Ptt9/SqlUr4uLi6Nq1K0lJSbzzzju89tprxMTEsHjx4tPmTk9Pp2vXrrRq1YqoqCi+/PLLIvyURURESp6hQ4cydOhQRo8ezdChQ52OI1IqLNt+iA170xgYV9fv78CZn18/euHF5S+y4dCGi3rMplWa8ljbx864zfr16xk1ahSLFy+mWrVqHDp0iIEDB9K3b1/uv/9+Jk6cyPDhw5kxY8YZj7Nq1SqWLVtG48aN6dChA4sWLWL48OGMHj2aWbNmERkZSUpKCqNGjeKHH36gfPnyPPvss4wePZqnn34agKpVq/Lzzz8D8Ntvv1G1alUeeeQRACZOnEiVKlXYv38/V199Nbfddhsul4vBgwczf/58qlWrRk5ODlWqVOH++++nQoUKPPLII6SlpTF16tRTZg4NDeWLL76gYsWKpKSk0L59e3r16nWuH7OIiIiIyEXz0ZIkwssFc1NMbaejFCu/bvacMm/ePPr06UO1atUAqFKlCkuWLOHDDz8EYMCAATz66KNnPU7btm2pXbs2AQEBxMTEkJSURMeOHQtss3TpUtatW0eHDh0AyM7O9r4GuOOOO057/Ndff50vvvgCl8vFrl272Lx5MwcOHKBz587Uq1ePtLQ0qlSpck7nbq3liSeeYP78+QQEBLBnzx727dtH+fLlz+k4IiIiIiIXw+9Hsvhu7T7u61iP0OBAp+MUK79u9s42AucLgoKCcLlcALhcLo4fP+5dFxIS4n0dGBhIbm7uH/a31tKtWzfvSFtaWhphYWHe9adrshISEvjhhx9YsmQJeXl53HjjjWRnZ19w7ilTpnDgwAESExMJDg4mMjKS7OxsNXsiIiIi4oiPl+3EZS1/al/X6SjFTtfsFYGrrrqKzz77jIMHDwJw6NAh4uLivNe4TZkyhU6dOgHu6+oSExMB+Oqrr8jJyTnr8cPCwkhLSwOgffv2LFq0iC1btgCQkZHBpk2bznqMI0eOULlyZcqVK8emTZtYunSp93jz589n+/bt3uwn1zxT7iNHjlCjRg2Cg4OZO3cuO3bsOGsWEREREZGikOOyTF2+k65Na3BplXJOxyl2avaKQLNmzXjyySe56qqraNGiBSNGjGDcuHFMmTKF6OhoJk+ezNixYwEYPHgw8+bNIy4ujiVLlhRqBGzIkCHceuutdOnSherVq/PBBx/Qt29foqOjueaaa9iw4ezXKXbv3p3c3FyaNWvGM888Q/v27QGoXr0648eP59ZbbyUuLs47DfTGG2/kiy++8N6g5UTuFi1aFMjdv39/Vq5cSVRUFB999BFNmzY9349RREREROSCrNibR0r6cQaWosct5OfX0zidNHDgQAYOHFhg2cyZMwtMsQSIiIhg6dKl3umXL774IgDx8fHEx8d7R9PeeOMN7z7Dhg1j0KBB3mNdffXVrFixAig4jTMpKalArSeeeMK7LiQkhNmzZ/9hH4Drr7+e66+/vsDyxo0bs3r16gLbnxgNBHjxxRdJS0ujWrVqLFmy5A+fR1pamvcuoSIiIrEc1R0AACAASURBVCIixeGHHTnUr16eDg2qOR3FERrZExERERERv7N4awrbjri4q31dAgJKz+MW8tPInoiIiIiI+I3cPBdvJ2zl9Z82UyXUcFtsHacjOcYvmz1rbal6WKL4D2ut0xFERERESqzN+9J4+LNfWb37CDe2qMV1VVMJCw12OpZj/G4aZ2hoKAcPHtQfzVLiWGs5ePAgeXl5TkcRERERKVFc1vLuvK30GLeQ3YezeLNfK8b1bUmFMqV7AMjvRvbq1KnD7t27OXDggGMZsrOzCQ0NLfTy89nHF2v4+/kVR43Q0FAyMjJOWVtERERE/mjHwQyeX5bNltQNXHt5BM/dEkX1sJCz71gK+F2zFxwcTL169RzNkJCQQMuWLQu9/Hz28cUa/n5+xVVDzyYUERERKZwt+9O5c/xSMrNdvHZHC26Oqa3LufLxu2ZPRERERET839YD6fSd4H4U2FPty3JLy9J7I5bTUbMnIiIiUozeffddADZu3EiTJk0cTiNSMm1PyaDfhKW4XJZpQ9qzZ32i05F8kpo9ERERkWI0ZMgQwD2lPz4+3tkwIiXQjoMZ9B2/lJw8y9TB7WkUEcae9U6n8k1+dzdOERERERHxTwcyXfQdv5Ts3Dym3NeOJjXDnI7k0zSyJyIiIiIiPm/34UxeWJ5Nrgni48HtaHZJRacj+TyN7ImIiIiIiE/buDeNO95dSlauZcp97biiViWnI5UIavZERERERMRnzd24n9veXkxOnotH24TSvLYavcLSNE4RERGRYpSY6L5r4MaNGwkLCyM2NtbhRCK+64NF23l25jqa1qzI+4Nas3HVMqcjlShq9kRERESKUevWrQu8t9Y6lETEd+W5LE9/uYaPluzgmmYRjL0zhvIhQWx0OlgJo2ZPRERERER8xtHsHF5LPMaagzsY2rk+j3ZvSmCAcTpWiaRmT0REREREfMKRzBx6v7OYbYfyePG2KO5oc5nTkUo03aBFRERERER8wr9nrWNbSgYjYkPV6F0EavZERERERMRxCRv3Mz1xN0M71+eKaoFOx/ELavZERERERMRRadk5PP6/32hYowLDuzZyOo7f0DV7IiIiIiLiqP/M3sC+o9lMfyCO0GCN6l0sGtkTERERERHHLN6SwsfLdnJvx3q0uqyy03H8ipo9ERERERFxRMaxXB7732oiq5ZjRLcmTsfxO5rGKSIiIiIijnj5u43sOpTFp0OvpGwZTd+82NTsiYiIiIhIsdt0OI8Plycx8Mq6tK1Xxek4fknTOEVEREREpFhl5+Qx8bdj1A4vy6Pdmzodx29pZE9ERERERIrVf75Zz95My5R+0ZQPUUtSVDSyJyIiIiIixebbNb/z4ZIdXBcZRIeG1ZyO49fURouIiIgUo8GDBwOQnJxMrVq1HE4jUrx2Hcrk0emraVGnEn0a5zgdx++p2RMREREpRuPHjwcgISGB+Ph4Z8OIFKOcPBfDp63CWhjXtxXbflvudCS/p2mcIiIiIiJS5F6Zs5FVO1N54bZoLqtazuk4pYKaPRERERERKVIJG/fz7rxt9Gt3GT2iL3E6TqmhZk9ERERERIrM4WwXIz79laY1w3i65+VOxylVdM2eiIiIiIgUiTyX5d3Vx8g6bnijX0tCgwOdjlSqqNkTEREREZGL4khmDpv2p7FpXxqb9qbx6+4jbDjk4uXe0TSsEeZ0vFJHzZ6IiIhIMTLGFHhvrXUoicjF89K3G/h4SSap387xLitfJpCGEWH0bhxM79g6DqYrvdTsiYiIiIjIeftlVypvJWzliqoBPHB1ExpHhNEoogK1w8tijCEhIeEP/5NDioeaPREREREROW/jftxMeLlgHmwZzPVXNXA6juSju3GKiIiIiMh5WbPnCD9u2M99HetRNkijd75GzZ6IiIiIiJyX13/cTMXQIO6Ki3Q6ipyCmj0RERERETln65KPMmfdPu7pWI+KocFOx5FTULMnIiIiIiLnbNxPmwkLCeLuuHpOR5HTULMnIiIiIiLnZOPeNGav2cugDpFUKqdRPV+lZk9ERERERM7JuJ82U75MIPd21KieL1OzJyIiIiIihbYn3cWs335nYFwk4eXKOB1HzkDNnoiIiIiIFNrXW49TNjiQ+zrVdzqKnEWRN3vGmEBjzCpjzEzP+3rGmGXGmC3GmE+MMWU8y0M877d41kcWdTYRERERESm8rQfSWfZ7HgPa16VKeY3q+briGNl7CFif7/2LwGvW2obAYeBez/J7gcOe5a95thMRERERER/x5twtBAegUb0SokibPWNMHaAH8J7nvQGuBqZ7NvkQuNnz+ibPezzru3q2FxERERERh209kM6XvyTT5dIgqoeFOB1HCsFYa4vu4MZMB/4DhAGPAIOApZ7RO4wxlwKzrbXNjTFrgO7W2t2edVuBdtbalJOOOQQYAhARERE7bdq0Ist/vtLT06lQoUKhl5/PPr5Yw9/PTzVKRm3VKBm1VcO3avj7+flajSFDhgDgcrkICAhg/PjxxVZbNfyjhlO131iVzZqUPJ6JtVxSpeSeR1Ecy0ldunRJtNa2PuVKa22R/AA9gbc8r+OBmUA1YEu+bS4F1nherwHq5Fu3Fah2phqxsbHWF82dO/eclp/PPr5Yw9/PTzVKRm3VKBm1VcO3avj7+flLDX8/P9Xw7dqJOw7Zuo/NtK99v7FEn0dRHctJwEp7mn6pKKdxdgB6GWOSgGm4p2+OBcKNMUGebeoAezyv93iaPzzrKwEHizCfiIiIiIichbWWF2dvoFqFMrpWr4QpsmbPWvu4tbaOtTYSuBP4yVrbH5gL9PZsNhD40vP6K897POt/8nSqIiIiIiLikIRNB1i2/RDDuzaiQkjQ2XcQn+HEc/YeA0YYY7YAVYH3PcvfB6p6lo8A/uFANhERERER8chzuUf16lYtx51tLnM6jpyjYmnNrbUJQILn9Tag7Sm2yQb6FEceERERERE5uy9/2cOGvWm83rclZYKcGCeSC6FvTERERERE/uBYbh6vztlE89oV6Rl1idNx5Dxo0q2IiIhIMYqNjQUgLS2NsLAwEhMTHU4kcmr/XbqTPalZvHhbNAEBevx1SaRmT0RERKQY/fzzz05HEDmrzBzLG/M306lRNTo2quZ0HDlPmsYpIiIiIiIFzN6ew+HMHB7r3tTpKHIB1OyJiIiIiIhXcmoW3+3IoVeLWjSvXcnpOHIB1OyJiIiIiAgAB9KOMeD9ZRjg4WsbOx1HLpCaPRERERERISX9GP0mLCU5NZsRsaHUrVre6UhygdTsiYiIiIiUcocyjvOn95ax63AmEwe1oUmVQKcjyUWgZk9EREREpBQ7nHGc/u8tY3tKBu8PbMOVDao6HUkuEj16QURERESklEo/bvnT+8vYeiCd9+5qTYeGesyCP1GzJyIiIiJSCh3KOM4rK7NJzoDxd8XSuXF1pyPJRaZmT0RERESklEjLzuGH9fuYtfp35m9KIc/lYsLA1sQ3qeF0NCkCavZERERERPxY5vFclibn8vFHK0nYdIDjuS4uqRTKXVfWpa7dy9VNI5yOKEVEzZ6IiIiIiJ+y1tJvwjJ+2XWMiIqp9G93GT2jL6HlpZUJCDAkJOx3OqIUITV7IiIiIiJ+asHmFH7ZlcqdTcrw/MCuBAQYpyNJMVKzJyIiIiLipyYs2EaNsBC61g1Uo1cK6Tl7IiIiIsXIWou1lrlz52KtdTqO+LG1yUdYsDmFQR0iCVajVyqp2RMRERER8UPvLdhO+TKB9G9X1+koJd6KvStIOpbkdIxzpmZPRERERMTPJKdm8fWvydzR5jIqlQ12Ok6JdvT4Uf6x4B9MPTgVl3U5HeecqNkTEREREfEzExduxwL3dIx0OkqJ9/KKlzmYdZB+VfsRYEpW+1Sy0oqIiIiIyBkdycph6vKd9Iy+hDqVyzkdp0Sbv3s+M7bM4O7md1M3pORNh1WzJyIiIiLiR6Yu30nG8TwGd6rvdJQS7cixI4xcPJKG4Q15oMUDTsc5L3r0goiIiEgxGjJkCADJycl8/PHHjB8/3uFE4k9yXZZJS7fToWFVmteu5HScEu2lFS9xMPsgr1/9OmUCyzgd57yo2RMREREpRhMmTCjwXs2eXExLf89l39HjvHhbtNNRSrT5u+fz1davGBw1mCuqXeF0nPOmZk9ERERExA9Ya/l2ew5NIsK4qnF1p+OUWJl5mbyy+BUahjfk/hb3Ox3nguiaPRERERERPzBv0wF2p1sGd66PMXqI+vn6/PDnHMo+xKiOo0rs9M0T1OyJiIiIiJRw1lrGz99GeIihV4taTscpsRJ2JbA8Yzn3Rt3LFVVL7vTNEzSNU0RERESkBLPWMvr7TSzeepA7m5ShTJDGc05nf+Z+Vu1fxdqDa1mbspbPfvyM1OxUDh87TOqxVNKOp1EruBb3R5fs6ZsnqNkTERERESnBxv64mXE/beHONpdybZWDTsfxKcnHk/lkwyesOrCKX/b/wp70PQAEBwRTwVSgZmZNwkPCqR1Wm8ohlQkPDSfiQATBgcEOJ7841OyJiIiIiJRQ437czJgfNtM7tg7P3xLF/PnznI7kM2Zum8l/fv8P/A7VylajZY2W9G3al5Y1WtKsSjMWLVhEfHz8H/ZLSEgo9qxFRc2eiIiIiEgJ9ObcLbz6/SZubVmbF2+LJiBAN2XJb+bWmVQLqsZHvT6iToU6pfKmNZrQKyIiIiJSwrw7bysvf7eRm2Jq8XKfFgSq0SsgIyeD5XuXE10umkvDLi2VjR5oZE9EREREpET5LimHqRs2cGOLWryqRu+UliYvJceVQ/OyzZ2O4iiN7ImIiIiIlBBf/rKHqRuO0yPqEl67vQVBgfpz/lTm7Z5HWJkw6ofUdzqKo/TbISIiIiJSAmxPyeCJ//1Gw/AAxtwZo0bvNFzWxfzd8+lYqyOBJtDpOI7Sb4iIiIiIiI87lpvHsKk/ExQYwAMtQghWo3daa1PWcjD7IJ0v7ex0FMfpt0RERERExMf955sNrNlzlFf6tKBqWf0Jfybzds8jwATQsVZHp6M4TjdoERERESlGK1eu9P6zdevWDqeRkuC7tXv5YHESd3eIpNvlESTsX+90JJ82f/d8YqrHEB4a7nQUx6nZExERESlGsbGxAKSlpXlfi5zOntQsHp2+mqjalfjH9U2djuPz9mXsY/2h9fy11V+djuITNAYsIiIiIuKDcl2W4VNXkeeyvNGvJSFBpftmI4Uxf898AOIvjXc2iI/QyJ6IiIiIiA/6YnMOiTsyeb1vS+pWLe90nBJh3q551K5Qm/qVSvcjF07QyJ6IiIiIiI9ZsPkAs7bn0LftZfRqUcvpOCVCVm4WS39fylV1rsIYPWge1OyJiIiIiPiUtOwcHp2+mlrlDc/ceLnTcUqMFXtXcCzvGFfVucrpKD5D0zhFREREitH48eMB2LhxI5s2bWLIkCEOJxJf85/ZG9h3NJsn24USGqzr9Apr3q55lAsqR+uausvtCWr2RERERIrR0KFDC7xXsyf5Ld6SwsfLdjKkc30alNvndJwSw1rLvN3ziKsVR5nAMk7H8RmaxikiIiIi4gMyjuXy2P9WU69aeUZ0a+x0nBJlT84e9mXuo3Odzk5H8Ska2RMRERER8QEvf7eR3Yez+GTIlZq+eY7WZK4BoFOdTg4n8S0a2RMRERERcdiKpEN8uCSJgVdG0rZeFafjlDhrs9YSVS2KamWrOR3Fp6jZExERERFx0PE8y6PTV1M7vCx/v66J03FKnJSsFHYc36EpnKegaZwiIiIiIg76YksO21NymHJfO8qH6M/zc7Vg9wIsVo9cOAX9NomIiIiIOGTVzsN8uz2Hvm0vpUNDTUE8FxsObWD6punM2jaLyoGVaVqlqdORfI6aPRERERERB+w+nMmDH68iPMTw+A3NnI5TImS7spm+aTqfb/qcNQfXUCagDNdGXktUVhTGGKfj+Rw1eyIiIiIixSw5NYu+E5ZyNDuHh1uFUDE02OlIPslay/Yj21m1fxWJ+xKZs3sOx3Ydo2F4Q/7R9h/0rN+TSiGVSEhIcDqqT1KzJyIiIiJSjPYeyabvhKWkZuTw3/vacXjrL05H8ilrU9Yy58gcpv84nV8O/MKRY0cACA8JJ6ZcDA92fpAW1VtoJK8Q1OyJiIiIiBSTfUfdjd7B9ONMvrctLS4NJ2Gr06l8xycbPmHUslEARLoiufrSq2lZoyUxNWKIrBjJvHnziKkR43DKkkPNnoiIiIhIMUg95qLvhKXsP5rNR/e2peVllZ2O5FNmbZvFc8ueI75OPNdyLTd2vdHpSCWenrMnIiIiIlLEDqQd46Xl2ew9ks0H97Qltq4enJ7fvF3zeHLhk7Su2ZpX4l8hLDDM6Uh+QSN7IiIiIiJF6GD6Mfq/t5SUbMvke9vQJlKNXn4r9q7g4XkP07RKU8ZdPY6QwBCnI/kNNXsiIiIixejdd98FYOPGjTRp0sThNFLUDmUcp/97y9h5KJO/tQqlXf2qTkfyKWtT1jLsp2HUqVCHt695m/LB5Z2O5FfU7ImIiIgUoyFDhgCQkJBAfHy8s2GkSKVmuhu97SkZTBzUhpzda5yO5FP25uzlnz/8k/CQcN7t9i6VQ3UN48Wma/ZERERERC6yI5k59H9vGVsPpDPhrtZ0aFjN6Ug+JTk9mTf3vUmgCWR8t/FElI9wOpJf0sieiIiIiMhFlJFjGTBxGZv3pfPuXbF0blzd6Ug+5XD2YYZ+P5Rj9hiTu03msoqXOR3Jb6nZExERERG5SI5m5/Dqymx2pWfxzp9i6dKkhtORfEpmTiZ/+fEv/J7xOw9Uf4AmVXTdalFSsyciIiIichGkpB9jyEcr2XHUxTsDWtO1maYm5pfjymHEvBGsPbiWMfFjMNuM05H8npo9ERERkWKUmJgIuO/GGRYWRmxsrMOJ5GL45rffeWrGGtKP5fLnmBC6Xa5GLz+XdfH0oqdZtGcRI+NG0uWyLiRsS3A6lt9TsyciIiJSjFq3bl3gvbXWoSRyMRzOOM7TX63l61+Tia5TiVf7tGDP+kSnY/mc0StHM3PbTIa3HM6tjW51Ok6poWZPREREROQ8/LBuH49/8Rupmcd5uFtjHohvQFBgAHvWO53Mt/x45Edm7JhB/2b9uS/qPqfjlCpq9kREREREzkFadg4TVh9jUfJKmtYM48O723J5rYpOx/I5BzIPMHXDVGakzqB7ZHcebfMoxug6veKkZk9EREREpJBy8lwMnZzI0t9zGXZ1Q4Zd3YgyQXp09Qku62Lp70v5bONnJOxKINfm0rJcS57r+BwBRp9TcVOzJyIiIiJSCNZa/jljDYu3HmRwVBkevlaPDTjhcPZhfjjyAy/+70V2p++mckhlBlw+gN6Ne7Pt522UCSzjdMRSSc2eiIiIiEghvLdgO9NW7OLBLg1pHfK703F8xrxd83hi4RMcPX6U2IhYhrUcxjV1r/E2eNvY5nDC0kvNnoiIiIjIWXy/bh/Pz17PDVE1GdGtMfPnq9nLceUwbtU4Jq2ZRNMqTflLyF/od20/p2NJPpo4KyIiIiJyBmuTj/DQtFVE167Eq31iCAjQTUb2Zezjvu/uY9KaSdze+Hb+e8N/qVWmltOx5CQa2RMREREROY3D2S7+8cFKwssGM+Gu1pQtE+h0JMetz1rP018/TXZeNi90eoEe9Xs4HUlOQ82eiIiIiMgpZB3P4/Wfj3E02zD9/jhqVAx1OpKjcl25vPPrO4zfP54G4Q14Nf5V6leq73QsOQM1eyIiIiIipzDy67UkHXUx4a7Wpf45eruO7uLxhY/z64FfaVe+HeN6jKNsUFmnY8lZqNkTERERETlJcmoW0xN3c03dIK65PMLpOI6x1jJjywxeWP4CgSaQlzq/RNkdZdXolRC6QYuIiIiIyEkmLdqOBa6LDHY6imNSs1N5eN7DPL34aa6odgWf9/qc6+td73QsOQca2RMRERERyedodg5Tl++iR9QlVCt7xOk4jtiQtYF/f/VvDh07xN9i/8bAywcSGKCb05Q0avZEREREitHgwYMBSE5OplYt3areF01dtpP0Y7kM6VyflM2rnI5TbLJys/hhxw/M2DKD5fuXU79Sfd7o+gbNqjZzOpqcJzV7IiIiIsVo/PjxACQkJBAfH+9sGPmD47kuJi1KIq5BVZrXrkTCZqcTFS1rLb/s/4UZW2bwbdK3ZORkUKdCHXqG9+TpHk/r2rwSTs2eiIiIiIjH178ms/doNi/cFuV0lCK3aM8inkt+jn0791E2qCzd6nbj5oY3ExsRy/x589Xo+QE1eyIiIiIiuEe5JizYRpOIMK5qXN3pOEVq7cG1/HXuXwkPCOfZuGe5NvJaygeXdzqWXGS6G6eIiIiICDBv0wE27E1jcOf6GGOcjlNk9mbsZdiPw6gSWoVhEcO4pdEtavT8lJo9ERERERFgwoJtRFQMoVcL/71xTmZOJg/++CCZuZm80fUNKgaW7ofF+zs1eyIiIiJS6iUdyWPRloPc3aEeZYL880/kPFcej85/lM2pm3nlqldoVLmR05GkiOmaPREREZFidPL0QGutQ0kkv2+TcihfJpC+bS9zOkqReTXxVebtnseT7Z6kY+2OTseRYqBmT0RERERKtd2HM1m+N497OtSjUtlgp+MUiQVpC/h0x6f8qdmfuLPpnU7HkWLin2PUIiIiIiKFNHFhEgB3d6znbJAisnDPQqYfmk7nOp15pPUjTseRYqSRPREREREptX7ZlcqUZTtoVzOQ2uH+9Vy5HFcOE1ZPYPzq8dQMrslLnV8iMCDQ6VhSjNTsiYiIiEiptCc1i/s+XEn1sBD6NvWvRy1sS93G4wsfZ93BdfSs35NOxzvp8QqlkKZxioiIiEipk34sl3s/WMGxnDwmDmpDxRD/aPZc1sVHaz+iz9d9SE5PZnT8aP7T6T+UCyzndDRxgEb2RERERKRUyXNZHpq6ik370ph0d1saR4SRvN7pVBduT/oe3tj3Bpt3bia+TjzPxD1DtbLVnI4lDlKzJyIiIiKlyn++Wc+PG/bz7E1XcFXj6k7HuSj2Zuzl9q9v53jOcZ6Ne5abG978h8d8SOmjZk9ERERESo2EXTl8sHY7g+IiuevKSKfjXBTWWp5d8iw5rhweueQRbml0i9ORxEfomj0RERERKRUWbk7ho3XHiW9Snad6NHM6zkUzc9tMFuxZwEOtHiIiOMLpOOJDiqzZM8aEGmOWG2N+NcasNcaM9CyvZ4xZZozZYoz5xBhTxrM8xPN+i2d9ZFFlExEREZHSZV3yUf48JZFa5Q3j+rYkKNA/xjxSslJ4ccWLxFSPoW/Tvk7HER9TlL/lx4CrrbUtgBiguzGmPfAi8Jq1tiFwGLjXs/29wGHP8tc824mIiIiIXJANe4/S/72llA8J4q+xoYSFBjsd6aJ5ftnzZOVkMbLDSAKMfzSwcvEU2W+EdUv3vA32/FjgamC6Z/mHwM2e1zd53uNZ39XoqlIRERERuQCb9qXRf8IyQoICmTq4PdXK+k9D9P2O7/l+x/c8EPMA9SvVdzqO+CBjrS26gxsTCCQCDYE3gZeBpZ7RO4wxlwKzrbXNjTFrgO7W2t2edVuBdtbalJOOOQQYAhARERE7bdq0Ist/vtLT06lQoUKhl5/PPr5Yw9/PTzVKRm3VKBm1VcO3avj7+flajS5duhRYP3fu3GKrXdpqJKe7eGF5FgHG8I+2odQsH1Aiz+NUy01Zw3PJzxEeFM7DNR8m0AQWS23VOPWxnNSlS5dEa23rU6601hb5DxAOzAU6AlvyLb8UWON5vQaok2/dVqDamY4bGxtrfdHcuXPPafn57OOLNfz9/FSjZNRWjZJRWzV8q4a/n5+v1WjVqpVt1aqVbdSokW3VqlWx1i5NNbbsT7OtR31vY//9vd28L61IapxOcfw7+8SCJ2zMhzF2w8ENxVpbNXwPsNKepl8qlkcvWGtTjTFzgSuBcGNMkLU2F6gD7PFstsfT/O02xgQBlYCDxZFPREREpLgkJiYCkJCQQHx8vLNh/NTeDBePjl+KtZapg9vTsIbvjcZciLVZa/lq/1cMiR5CkypNnI4jPqzImj1jTHUgx9PolQW64b7pylygNzANGAh86dnlK8/7JZ71P3k6VRERERGRQtl1KJMXl2djgoKZOrg9jSLCnI50RrmuXA7nHubo8aOUCypHUEDBP89d1sXh7MPsz9zP/sz97MvcxycHP6FBpQYMjR7qUGopKYpyZO8S4EPPdXsBwKfW2pnGmHXANGPMKGAV8L5n+/eBycaYLcAh4M4izCYiIiIifuZIVg6DJi3nuMsy/b52NKnp241eniuPe767h1X7V/H01KcBKBNQhvLB5SkXXI6srCzS/ptGjiunwH5lTVme7fAsZQLLOBFbSpAia/astauBlqdYvg1oe4rl2UCfosojIiIiIv4rJ8/Fgx//zM5DmTwcG0qzSyo6Hemspm2cxqr9q7im4jW0bNySzNxM909OJlm5WSTvTSaqXhQR5SK8PzXK1fg/9u47PMoq7//4+6SHJCQQIFTpvROUZgkWUFDsDfujxPqs+1t3V1139dld3VVcdS27KqhYAUVdccECAkOvoQRCC6GTBoSEhJBkMnN+fyQiiIQJzGSSyed1XXNN7no+9x8Qvpxzn8OGFRvo07SPv+NLHVAj7+yJiIiIiPiKtZb/+zqNhekHGH9DH5oVZfg70mnlFufy+prXGdZyGGNCxjC85/CTznE4HCQlJp20f5PZVAMJJRAEzkIjIiIiIlIvTVq8k0+W7+aBizpy08A2/o7jkfErx+N0OXlq0FNoaWnxFfXsiYiIiNSgxMREAAoLC4mJiTk2O6ecmbmbc3h25kZG9kzg9yPrxsyUi/Yt4vud3/NIv0do07ANGdT+nkipm1TsiYiIiNSg1atX+ztCwNhT6Ob5uWvo0bIhr9zcj6Cg2t9DVlJewnPLnqNdapPclQAAIABJREFUw3bc0+sef8eRAKdiT0RERETqnNzCEv6ZUkJ0RBjv3HkuDcLqxj9rJ6ROYG/RXt4d8a5m0xSfqxt/KkREREREgOKycj5auou3F2ynyGn5cty5NI+N8Hcsj2Q7s5mUNomrOlzFeS1OmpxexOtU7ImIiIhIrVfidPHxsl28NT+DA0VlXNilKRfHF9KrVay/o3nEWsunBz8lMiSSxwY+5u84Uk9Uq9gzxjQC2lSuoSciIiIi4lMlThezdzr53eJ57C8sZVineN66tAsD2zXG4XD4O57H/rv9v2wr3cbTQ54mPjLe33GknjhtsWeMcQBjKs9NAXKNMYuttb/xcTYRERERqceWZBzg95+nsvdQGYPaN+aNW/szqEPdKpRyi3P5eNPHTN08lfbh7bm+8/X+jiT1iCc9e7HW2sPGmPuAD621zxhj1LMnIiIiIj5RXFbOC99u5oOlu2jfJIrfnxvBQ9cP8Xesatl2aBvvp73PzB0zcVs3l7W9jKHOoQQZLXMtNceTYi/EGNMCuAl4ysd5RERERKQeW7Uzj99OW8fOg8XcPbQdj1/ejeVLFvo7lsdWZq/kzdw32fj1RiJDIrmpy03c3uN22sS0qVPDTiUweFLs/QX4HlhsrV1pjOkApPs2loiIiIjUJyVOF1M3l/L990tp3SiSqcmDGVyHhmxuL9jO+JXjWbxvMTFBMTzS7xFu7nozcRFx/o4m9dhpiz1r7TRg2nHb2wENNhYRERERrzhYVMotE5aRnlvObYPO4Q+juhMVXjcmjT9cdpg3177J1M1TiQiJ4LcDf0uLnBaM6DvC39FEPJqgpQvwJpBgre1ljOkDjLHWPuvzdCIiIiIS0Nxuy2PT1rErr5jfJIbzq2t7+zuSR9zWzbSt03h99evkl+ZzfZfreaTfI8RHxuPY7/B3PBHAs2GcE4HfAW8DWGtTjTGTARV7IiIiInJW3lm0HceW/fz16p60Kd3p7zge2XpoK+OzxrNv9z4SExJ5/NzH6R7f3d+xRE7iSbHXwFq7whhz/L5yH+URERERkXpize5DjP9uC1f0as7tg9syf/5Of0c6rbSDaSTPSgY3/OOifzCi7Qh+9u9kkVrDk2LvgDGmI2ABjDE3AFk+TSUiIiIiAa3gqJP/nbKG5rERPH99nzpRMK3fv577Z99Pw/CG3NfwPka2G+nvSCJV8mShj4epGMLZzRizD/g18KBPU4mIiIgEKGst1lrmzZuHtdbfcfzCWssTX6SSXVDC67f2JzYy1N+RTmtt7lrGzR5HbHgsk0ZOokloE39HEjktT2bj3A5caoyJAoKstYW+jyUiIiIigerjZbv4dkM2T17Rjf7nNPJ3nNNKyUnhoR8eommDprwz4h2aRzVnC1v8HUvktE7bs2eMedQY0xAoBl4xxqw2xmguWRERERGptl2HXfx15iaSujZl3AUd/B3ntLaWbOXBHx4kISqBSSMn0Tyqub8jiXjMk2Gc/2OtPQyMAOKBO4DnfZpKRERERALO4RInb64tpVGDUF66sS9BQbX3Pb28kjy+2PoFb+W+RavoVrw38j2aNmjq71gi1eLJBC0//ikcBXxorU0zdeENWhERERGpNZZsO8DvPk8l96jlk/v6Ex8d7u9IJ3BZF6uyV7EkcwmLMxez6eAmLJbWYa15d+S7NI5o7O+IItXmSbGXYoyZBbQHnjTGxABu38YSERERkUBQXFbO899u5sOlu2jfJIonz4tgSMd4f8c6ptxdzrPLnmXmnpmU7C4h2ATTp2kfHur3EMNaDmP/hv0q9KTO8qTYuxfoB2y31hYbY+KBe3wbS0RERCQwJScnA5CZmcnkyZOZMGGCnxP5zqqdeTw2bR27DhZz99B2PH55N5YvWejvWCd4YcULfJH+BedFncet597KoBaDiAmLOXbcYRz+CydyljyZjdNtjNkBdDHGRNRAJhEREZGANXHixBO2A7HYK3G6mLq5lO+/X0rrRpFMTR7M4A61pzfvR1M2T2Hqlqnc3fNuEosSSWqb5O9IIl512mLPGHMf8CjQGlgLDAaWAhf7NpqIiIiI1DVut+WhT1Yzd2c5tw06hz+M6k5UuCeDyWrWkn1LeGHFCyS1TuLXA37NwgW1q8dRxBs8mY3zUeBcYJe1djjQH8j3aSoRERERqZPeXbSDuZtzua1bGM9d27tWFnrb87fz2/m/pWNcR56/8HmCg4L9HUnEJzwp9kqstSUAxphwa+1moKtvY4mIiIhIXbN2Tz4vfLeZkT0TuLRt7SvyAI64jvDI3EcIDQ7l9YtfJyo0yt+RRHzGk2JvrzEmDvgKmG2MmQ7s8m0sEREREalLCo46eWTyahIaRjD++r7UxpW6nC4n7+x/h5wjObx28Wu0jG7p70giPuXJBC3XVv74f8aYeUAs8J1PU4mIiIhInWGt5ckvU8kuKOGzB4YQ2yDU35GAilzF5cXkHc0jrzSPTzd/yrbSbTx/wfP0bdrX3/FEfO6UxZ4x5lygibX22x/3WWvnG2NGAb2BlBrIJyIiIiK13MfLd/PN+myeuKIbA85p5JcMbusmdX8qP+z6gXlZ8/jb538jrySPUlfpCeddHns5ozuM9ktGkZpWVc/eC/zyenppwCQ0G6eIiIhIvbf7sItnf9jIRV2aknxBhxpt22VdLM9azuxds5m7ey77j+4nJCiE9qHt6de8H/ER8TSKaETjiMY0jmhMQlQC+9buq9GMIv5UVbEXY6096d08a+0uY0wTH2YSERERkTrgSGk5/15XSlxkGC/f1JegIO+/p7c9fzspR1LI2pxFfkk++aX5HCo9RH5JPqk5qRzZfYSI4AjOb3U+l7a9lAtbX0jKkhSSzk/6xftlmkyvZxSpraoq9qrqg2/g7SAiIiIiUrs5XW52HjjC1pwituYUsjB9PzlHLJPH9Sc+OtyrbRWWFfLa6tf4dMunWCwcqNgfExZDXHgcjcIb0SOyB2PPHcuwVsOIDIn0avsigaCqYu8HY8xzwB+ttRbAVEyr9Gdgbk2EExEREZGaV+5ysyuvmPScwmOF3ZrtxeTO/g6nywJgDLRt3IA7eoQxpGO819q21vLdzu8Yv3I8B48e5JZut9CuoB0jLhhBbHgsoUE/Tf7icDhIapvktbZFAk1Vxd5jwDvANmPM2sp9fYFVwH2+DiYiIiIiNWdrTiETU0t5Yd1CMvYXUVbuPnasTeNImkQGcVViO7o2j6Zzsxg6NYsmIjQYh8PhtQy5zlzun30/S7OW0iO+B29c/AY9m/TE4XDQJFJvEYlU1ymLPWvtEeBWY0wHoGfl7jRr7fYaSSYiIiIiPudyWyYs2M4rs7cSYtyc2yGcCzo3oXOzaLokVBR1UeEhFb1oSd3Ovj3rIq8kj4LSAgpKCzhcdpiC0gLS89P5KPMjIkIj+MOgP3BTl5sIDgr2whOK1F+erLO3HVCBJyIiIhJgMvYX8dtp61izO5/LezZnVLPDjBl5nlfbOFp+lDW5a1ietZzlWcvZeHAjdrf9xXMHNBjAP0b/g6YNmno1g0h9ddpiT0RERES8Z9WqVce+Bw4c6JcMbmt5d9EOxn+3mYjQYF69pR9j+rZk/vz5Ht+joLSA73d+z4pDK9i0bhNhQWGEBoUSFhxGWHAYK/JX8MF3H7Bu/zqcbichJoQ+TftwacNLSeyWSGx4LHHhccSGxVb8HBHH6iWrVeiJeJGKPREREZEalJiYCEBhYeGxn72trNzNkowDzN3tZO+yE1fSssDHK0rYcmgjF3drxt+v601CwwiP7mutZVXOKr5I/4LZO2dT5i7DYLBrT+6pMxi6BXfjtu63MajFIAY0G0CD0AYVw0G7J3nhKUXkdE5Z7BljGld1obU2z/txRERERORMOF1ulmYcZEZqJt+n5VBw1FlxYOOGk86NDIHxN/ThxsTWVEy2XrX9xfuZVTCLF//zIrsLdxMTGsO1na/lus7Xkbs+lwsuvACn20mZuwyny4nT7WT1stWMumSUtx9TRKqhqp69FCr+88cA5wCHKn+OA3YD7X2eTkRERESqlLIrj0kbSvl/C37gULGT6PAQRvRI4Mq+LTi8M42hw4aedM3aFUsZMbBNlfc94jzCnN1zmLl9JsuyluG2bgYmDOSBvg9wadtLj61rl0suwUHBBAcFE8FPPYQNgrUss4i/VTUbZ3sAY8xE4D/W2m8qt68ArqmZeCIiIiJyKtNW7eH3X6QSHgQje7dkdO8WXNilKRGhFbNYOrI30Szm5CGaYcG/3JtXbstx7HEwc/tMHHsclLhKaBXdint73Uvzg8256bKbfPo8IuJdnryzN9haO+7HDWvtt8aY8T7MJCIiIiKn8eXqvfz+i1TO79SE29sWM/LS/qe9xulyMn/vfL7O+5p5S+ZRVFbEEecRipwV3/sK9nF091HiwuO4utPVjO4wmn5N+2GM8ep6eiJSMzwp9jKNMX8EPq7cvg3I9F0kERERkcA1YcIEALZs2cLWrVtJTk6u9j2mr93Hb6etY2jHeCbeOZBlixdWef7WQ1v5attXzMiYwaHSQ4SaUBrtbURUWBTRodFEhUbRNLIpCeUJ3HrerQxtNZTQoNAzej4RqT08KfZuBZ4B/kPFO3wLKveJiIiISDXdf//9J2xXt9j7el0m/+/TtQxqH887d557bMjmj8pcZRW9dWVFLCxcyFsz3iLtYBohQSEMbzOcaztdS1l6GZcMv+SkezscDi5qc1H1H0pEaiVPFlXPAx41xkRZa4/UQCYRERER+QUzU7P4f5+uZWC7xrx9Zz8c+2bx6ZZP2XFgB66pLo44j1DuLj/hms6NOvP4uY8zusNoGkU0AsCxzeGH9CJS005b7BljhgLvANHAOcaYvsD91tqHfB1ORERERCqsyi7nzVlr6HtOBCOHbOWGGc+QdSSLdg3b0TGiIx1adyAqtGJYZoPQBkSFRnE44zB3jLjDo+UVRCTweDKM8xVgJPA1gLV2nTHmQp+mEhERERGOlJYzZ3MuM9ZlMic9l1YdVpPVYDH/XFPEgGYDePK8J7mozUUsmL+ApMFJJ13v2OtQoSdSj3lS7GGt3fOzvyhcvokjIiIiUr+VuiwzU7OYuT6TuZtzKbX5xLWYT1THZRQEWS5rdRl39biL3k17+zuqiNRynhR7eyqHclpjTCjwKLDJt7FERERE6odl2w+yNaew8lPE2t3FlLlWE9/QSfeey9nlnIW1boZGDeapEU/RJqbqxdBFRH7kSbH3APAq0ArYB8wCHvZlKBEREZH64pYJywCIiQihS0IMQ1uV06JbOvOyPyejtJgrO1zJg/0eJCMlQ4WeiFSLJ7NxHqBibT0RERER8bIJd/UmKDyTnJIMNuYt4IcdP5Cyp5hLz7mUh/s9TKdGnQDIIMPPSUWkrvFkNs72wP8C7Y4/31o7xnexRERERALP3M05J+373cqrcFs3AI0jGtMlogtPXPwEPeN71nQ8EQkwngzj/Ap4F/gv4PZtHBEREZHAtCnrMP87ec1J++/vcz/dG3enR3wPmjVoxvz581XoiYhXeFLslVhrX/N5EhEREZEAlVtYwr3vryQmIvSkYw/109LFIuIbQR6c86ox5hljzBBjzIAfPz5PJiIiIhIAylyWcR+mcKjYyTt3DfR3HBGpRzzp2esN3AFczE/DOG3ltoiIiEj9VXIYdi6E3UvpvDMdCj4HlxPcTnCVYV1Otu1pR2r+hbx1eyK9WsX6O7GI1COeFHs3Ah2stWW+DiMiIiJSq1kX7E2BjLkVn70rwF0OweE0DYqAwigIDoHgMFwmlNwjLnIPN+OJy7sxsmdzf6cXkXrGk2JvAxAH5Po4i4iIiEjtlfoZwxb/BuYXVmy36AdDfwUdL4Y257Fk0VL6DxrG7I05zEjNZFH6AcrdlqTWITxzYYdjt3n77bcB2LJlC127dvXHk4hIPeFJsRcHbDbGrARKf9yppRdERESk3lg7Bb56kOKG3Yi97LfQIQmimgBQVu7m2w1ZTEopYePsHyhzuWkVF8m957fnyj4tOZC+GmPMsVslJycD4HA4SEpKqvlnEZF6w5Ni7xmfpxARERGprdZOhq8egg4Xsa7VQ1zYeyQATpebL1fv5bU529iXf5TGEYY7h7RjdJ8W9GsTd6zAc2wzVd1dRMRnqiz2jDHBwNvW2m41lEdERESk9jiu0OPWqbgXL6fc5eartZm8Nied3XnF9Gkdy7PX9MJmpXHx8B7+TiwickyVxZ611mWM2WKMOcdau7umQomIiIj43ZpPYPrDxwo9V3AESzLL+fMrC9hx4Ag9Wzbk3bsGcnG3ZhhjcGRv9HdiEZETeDKMsxGQZoxZARz5cafe2RMREZFA1TxrDjher3g379Yp7Chw87tpS1m1q5TuLcJ5+45ERvRIOOFdPBGR2saTYu9PPk8hIiIiUhu43bDibbpuqSj03DdP5sMV2Tz/3WbCgoMY1zuMJ289n6CgMy/yUlJSgIrZOGNiYkhMTPRSeBGRE5222LPWzq+JICIiIiJ+lZMG//017F1BXuOBHL38XX73QSpLtx8kqWtTnr+uD5vXLDurQg9g4MCBJ2xba8/qfiIip3LKYs8Ys8hae74xphA4/m8hA1hrbUOfpxMRERHxtbJiWDAelrwO4Q2x17zJ66lNmPbGSowxvHB9b24a2AZjDJv9nVVEpBqq6tm7E8BaG1NDWURERERqVKO81fDmo3BoJ0Xdb2Zao2Q+m3+UTVmHGdYpnheu70PrRg38HVNE5IxUVexNAxKNMXOstZfUVCARERERnykvg+z1sGc5bJ9H3/RZ5Ddoy/iGf2fymrZADv3PieOeXmH8aeygsx6yKSLiT1UVe0HGmD8AXYwxv/n5QWvty76LJSIiIuIFzhLYsYB2GVMp2fp3QnPWEOwqBSA7KIEpzut5M28M3Vo34ckrWjC6TwtaN2qAw+FQoScidV5Vxd4twDWV52gop4iIiNQN5aWQMRfS/oN700yCnEW0ssGk2Xascl9MirsLe6N706h5W5rZfGZfM4y28VH+Ti0i4nWnLPastVuAF4wxqdbab2swk4iIiEj1WEvjg6vgP1Nh80woPYwzLJavnecxL2go7ibduGBgX/onRHNjsxhiI0MBcDgcKvREJGB5ss7eXGPMWKDd8edba//iq1AiIiIi1TJ/PH3W/w3CY6H7VcwJHsrDSxvSrlkc79w1kG3rVpB03jn+TikiUqOCPDhnOnA1UA4cOe4jIiIi4n+5m2DBi+Q2HYbrt+n8Ofhh7l3ciKFdWvD5g0M1m6aI1Fue9Oy1ttZe7vMkIiIiItXldsH0RyA8hvUdkvn443XM27Kfe89vzx9GdSdYk6yISD3mSbG3xBjT21q73udpRERERKpjxQTYt4qDI9/gmTlhZBUf4NlrenH74Lb+TiYi4neeFHvnA3cbY3YApYABrLW2j0+TiYiIiFTl0C6Y8xcKWicx8ocEjpQ6+eCeQZzfuYm/k4mI1AqeFHtX+DyFiIiISHVYCzN+TbmFMTtvoEHDUH7TP0SFnojIcU47QYu1dhcQB1xV+Ymr3CciIiLiF3bdFMiYy5+P3kRC60589fAwWkZ7Mu+ciEj9cdqePWPMo8A44MvKXR8bYyZYa1/3aTIRERGRX3L0EMULHmejuwtH+9zNR9f3ITwk2N+pPDZu3DgAMjMzadmypZ/TiEgg82QY573AIGvtEQBjzAvAUkDFnoiIiNSo/YWl2JVvE+IqJn3Q33hxVD+MqVszbk6YMAGoWNA9KSnJv2FEJKB5UuwZwHXctqtyn4iIiEiN+W5DNrO/fJeX3EvZ2vNRxo6+zN+RRERqNU+KvUnAcmPMfyq3rwHe9V0kERERkZ8UFDt55usNbFi3gq8i/kVeZHu6XPdHf8cSEan1TlvsWWtfNsY4qFiCAeAea+0an6YSERERAeZtzuXxL1IxR3KZ1fCfRIVGs77XUwwJCfN3NBGRWs+TCVoGA2nW2tWV2w2NMYOstct9nk5ERETqpcMlTt5dX8rC71bSu1kYnzaaQIO8Q3DnN5SmH/Z3PBGROsGTOYrfBIqO2y6q3CciIiLiddtyC7nytUUs2lfOQxe156tWH9Igdy1cPxFaDfB3PBGROsOjCVqstfbHDWut2xjjyXUiIiIi1bJg634enrya8JBg/jAognGhn8Gmr2HEs9D9Kn/H84qfzx563D+zRES8ypOeve3GmF8ZY0IrP48C230dTEREROqXj5bt4p73V9IqLpLpjwzjgpJ5sOgVSLwHhjzi73giInWOJ8XeA8BQYB+wFxgEJPsylIiIiNQfLrfl/75O409fbSCpS1M+f3AorfKW02Xrv6HjJTDqH1DH1tITEakNPJmNMxe4pQayiIiISD1TWOLk1dWlpB7Yyb3nt+cPV3QjeOUEmP00xQ1aE33j+xCst0dERM6E/vYUERERvygscXLjW0vZetDFc9f24raeDWDqzZA+CzqPZF3T2xgW0dDfMUVE6ixPhnGKiIiIeN3z325mS04hvx4Qzm2Nt8KbQ2D7/Iphm2M/xRkW6++IIiJ12imLvcqJWDDGDKu5OCIiIlIfLNl2gE+W7yZ5SEuuy38fPrkBoppCsgPOG6d39EREvKCqnr17Kr9fr4kgIiIiUj8Ul5Xz+Jep9G1czuN7H6H1vv/CoAdg3DxI6OHveCIiAaOqd/Y2GWPSgZbGmNTj9hvAWmv7+DaaiIiIBKIXv99Cdl4hM9u+QdCBdNb3+iO9r/idv2OJiAScUxZ71tpbjTHNge+BMdW9sTGmDfAhkABYYIK19lVjTGPgU6AdsBO4yVp7yFSsMPoqMAooBu621q6ubrsiIiJSe63amcf7S3bwRctPaZizAq5/l4MHm/g7lohIQKpyghZrbba1ti+QBcRUfjKttbs8uHc58Ji1tgcwGHjYGNMDeAKYY63tDMyp3Aa4Auhc+UkG3jyD5xEREZFaqsTp4vefp/LbqO8ZkDcTLnoCet/g71giIgHrtLNxGmMuAtKBfwH/BrYaYy483XXW2qwfe+astYXAJqAVcDXwQeVpHwDXVP58NfChrbAMiDPGtKjm84iIiEgt9coPW+mU5+Ch8o+g53WQ9MTpLxIRkTNmrLVVn2BMCjDWWrulcrsLMMVam+hxI8a0AxYAvYDd1tq4yv0GOGStjTPGzACet9Yuqjw2B3jcWrvqZ/dKpqLnj4SEhMSpU6d6GqPGFBUVER0d7fH+M7mmNrYR6M+nNupG22qjbrStNmpXGzXRdlpWEd+k7uLL8D9TGnMOa/s9hzs4vM49hzfaGD58+AnH582bV2Ntq43AaCPQn6+uteFvw4cPT7HWDvzFg9baKj9Aqif7qrg+GkgBrqvczv/Z8UOV3zOA84/bPwcYWNW9ExMTbW00b968au0/k2tqYxuB/nxqo260rTbqRttqo3a14eu2S5zldsxfPrI5z7Szrpe6W3s42+ttVLW/trVBxVwGxz412bbaCIw2Av356lob/gassqeol6qajfNHq4wx7wAfV27fBqyq4vxjjDGhwBfAJ9baLyt35xhjWlhrsyqHaeZW7t8HtDnu8taV+0RERKQOe2fWGv7ifInGYSUEjf0aYhL8HcmvBgwYAEBhYSExMTF+TiMigcyTYu9B4GHgV5XbC6l4d69KlUM03wU2WWtfPu7Q18BdwPOV39OP2/+IMWYqMAgosNZmefIQIiIiUgs5S8ie8zq3LX+FhkHFBN04BZr39ncqv0tJSQHA4XCQlJTk3zAiEtBOW+xZa0uBlys/1TEMuANYb4xZW7nvD1QUeZ8ZY+4FdgE3VR77hoplF7ZRsfTCPYiIiEjd43bBuqnYec/R/PA+lgT1o7z3nVzY9Qp/JxMRqVc86dk7I7ZiohVzisOX/ML5looeRBEREamLrCX+wHJ48wnYv4ncmB78uuxu7hp7JxEHNvs7nYhIvXPapRdERERETitzDbw/mt4b/gZuJ5kj3uKCvD/RuNelXN6rub/TiYjUSyr2RERE5MwdzoL/PAgThsP+LWzt/ADl9y/hwdVtiI4I5S9jevo7oYhIvXVGwziNMcnW2gneDiMiIiJ1Q5CrFBwvwOJ/grschv0KLniMzGVrmLd0L+v2FvDarf2Jjw73d1QRkXrrTN/ZO9W7eCIiIhKojh6CnDTIXMN5K16B0oPQ4xq47M/QqB0AWUVuXlq2lct6JHBVnxb+zVtLJSYmAj8tvfDj7JwiIt52RsWetfZtbwcRERGRWsRayJhDux2fQuZbkLMBCvYcO1wW04mIsR9D26HH9rnclvc2lBIREsRz1/SiYhUm+bnVq1f7O4KI1BOnLfaMMa2B14HzAUvFOnuPWmv3+jibiIiI+MORgzDjUdj0X9oSBE06Q5tBcO69kNAbmvdi9apNJB1X6AG8NT+D9Hw3/7ixN80aRvgpvIiI/MiTnr1JwGTgxsrt2yv3XearUCIiIuIn6bNh+sNQnAeX/YWFJd248JKRJ59nTlxKYdLiHbz4/RbOax7M9QNa1VBYERGpiiezcTa11k6y1pZXft4Hmvo4l4iIiNSksiMw4zfwyQ3QIB6S58GwR3EHn36ClQ+X7uTP/93IyJ4JJPcJ1/BNEZFawpNi76Ax5nZjTHDl53bgoK+DiYiISM2IObwV3r4QVr0LQx6BcfOgeW+Prv1k+S6enp7GZT0SeP3WAYQEqdATEaktPBnG+T9UvLP3ChXv7C0B7vFlKBEREakh6z9nwOrHoWFLuPNr6HCRx5dOXbGbp/6zgUu6NeNfYwcQFqLle0VEapPTFnvW2l3AmBrIIiIiIjVp80z4MpmC2O7EPfANRMZ5fOnCvU7eS1tPUtem/Pt2FXoiIrXRKYs9Y8zTVVxnrbV/9UEeERERqQnb5sC0u6Flf9a3f4wLPCz0CoqdfLpqN+9tKOP8zk146/ZEwkOCfZtVRETOSFU9e0d+YV8UcC8QD6jYExERqYt2Loa6BdLqAAAgAElEQVSpt0GTrnD757iWr6vy9MMlTman5TBzfRYL0/fjdFl6NQlm4p0DiQhVoSciUludstiz1r7048/GmBjgUSre1ZsKvHSq60RERKQW25sCk2+GuDZwx38gstEvnlbucvNdWjbvri4hbfYPlLnctIqL5J5h7RnduwV529ao0BMRqeWqfGfPGNMY+A1wG/ABMMBae6gmgomIiIh3RRXthI+fgah4uHM6RJ+8kpLLbZmRmsmrP6Sz/cARGoUbbh/cjiv7tqB/m7hjyyo4MjTrpohIbVfVO3svAtcBE4De1tqiGkslIiIi3nUgnb7rnobI6IpZNxu2POGw2235ZkMW//whnW25RXRrHsNbtw8gbP9mLh7ew0+hRUTkbFTVs/cYUAr8EXjquAVSDRUTtDT0cTYRERGppvziMjLyXcTu/mkgTmhRJl1m3oDLwsbLPqK0sCEU/nR8WVY5f391IVtyCuncLJp/jR3AFb2aExRkcDi2+OMxRETEC6p6Z09zKIuIiNQh1lpue2c5aZklsGwJAHEUMi3sL5SYQ9xS9ic2fpINZJ90bYemIbx2a39G925BsBZGFxEJCJ4sqi4iIiJ1wML0A6RlHuaqjqFcd2E/gp1H6DvvLqIP7Sd1+Htcnh3J7/r0Oem6zRvWk3ztRSryaoi1FgCHw0FSUpJ/w4hIQFOxJyIiEiAmLtxO05hwxnQMZnjHOJhyP+Sth5s+on/3KylwOEjq2uyk60xWsAo9EZEApKGaIiIiASAts4CF6Qe4e2g7Qo2Frx6AjLlw1WvQ/Up/xxMRET9QsSciIhIA3lm4gwZhwdx+3jl02vYObPgCLv0zDLjD39FERMRPNIxTRESkjsvMP8p/12Vy55B2xK59i9h9M2HIIzDsUX9HExERP1LPnoiISB03afEOLPA/57eDlPc5FNcbLvsrGL2HJyJSn6lnT0REpA47XOJkyoo9jO7dgtaR5ZCXQX6722gUpP/Pra2Sk5MByMzMZPLkyUyYMMHPiUQkUKnYExERqcOmLN9NUWk5yRd2gJwNABTGdPBzKqnKxIkTT9hWsScivqL/9hMREamjyt2WSYt3MrRjPL1axULWOgCKojv6OZmIiNQGKvZERETqqOVZ5WQfLqno1YOKYi+6OWXhjfwbTEREagUVeyIiInWQtZZvdzjpmhDDRV2aVuzMSoUWffwbTEREag0VeyIiInXQgvQD7C2yjLuwA8YYcB6F/ZuhRV9/RxMRkVpCxZ6IiEgdY63l7fkZxIUbxvRtWbEzZyNYl4o9ERE5RsWeiIhIHfPOwh0syTjIqPahhIVU/irPWlvxrWJPREQqqdgTERGpQ2alZfO3bzcxuncLLm173ApK2akQEQexbfwXTkREahUVeyIiInXEhn0FPDp1LX1ax/HSTX0JMuang1nrKnr1jt8nIiL1moo9ERGROiDncAn3fbCKRg1CmXhnIhGhwT8ddDkhJ01DOEVE5AQhpz9FRERE/Km03HLfB6soLHHy+YNDaRYTceIJ+zeDq0zFnoiInEDFnoiISC3mdlsmrC8lLbeYiXcOpHuLhieflJVa8a1iT0REjqNhnCIiIrXYi7O2kJLj4qnRPbike8Ivn5S1DsKioXHHmg0nIiK1moo9ERGRWmpGaiZvOjJIahPC/wxrd+oTs9ZB894QpF/rIiLyEw3jFBERqYV2HTzCk1+sp/85cdzerQxzqlk23S7IXg/9b6/ZgHLGVq1adex74MCBfk4jIoFMxZ6IiEgtU1bu5n+nrMEYeP3W/mxbt+LUJ+dtB+cRva9XhyQmJgJQWFh47GcREV/QeA8REZFa5oXvNpO6t4DxN/ShdaMGVZ+cta7iW8WeiIj8jIo9ERGRWmTOphzeXbSDu4a05fJeLU5/QdZaCA6Hpl19H05EROoUFXsiIiK1RFbBUR6bto4eLRry5KjuHl60DhJ6QHCob8OJiEido2JPRESkFnC5Lb+asoaycjdvjO1PRGjw6S+ytmKNPQ3hFBGRX6AJWkRERGqBrzKcrNxZzD9v7keHptEeXRNRkgsl+Sr26pgJEyYAsGXLFrZu3UpycrKfE4lIoFKxJyIi4kdl5W4mL9/FjAwnNya25pr+rTy+Nrpoe8UPKvbqlPvvv/+EbRV7IuIrKvZERET8wOly83nKXt6Yu419+Ufp1jiIP1/ds1r3iCnMABMMzap3nYiI1A8q9kRERGpQudvy2co9vD4vnT15R+nXJo6/X9cb174NNAir3q/l6KLt0LQbhEb4KK2IiNRlKvZERERqyOyNOTy16Ci5xan0aR3LX8b0IqlrU4wxODJNte8XXbQdul/hg6QiIhIIVOyJiIj4mLWWfzsyePH7LbSONrxz50Au6d4MY6pf4B1TmE142SG9ryciIqekYk9ERMSHSstdPPnFer5cs49r+7fiiiaHuLRHwtnfOGtdxXeLPmd/LxERCUhaZ09ERMRH8o6Ucfs7y/lyzT4eu6wLL9/Ul7Dgs+jNO15WasV3897euZ+IiAQc9eyJiIj4QGaRm6f/tZicwyW8MbY/V/Zp6d0GstZSHNmSBuEx3r2viIgEDBV7IiIiXrYwfT9/XXaUqIhwPr1/CP3axHm3gfIyyFxDYUxHGnj3ziIiEkA0jFNERMSLPlq2i7snraRJZBDTHxnm/UKvYC+8PwoO7+Ng/LnevbeIiAQU9eyJiIh4QbnLzbMzN/H+kp1c3K0ZN7YuolVcpHcbyZgLX9xX0bN34wfk7o+jh3dbEBGRAKKePRERkbNUWOLkvg9X8f6Sndx7fnsm3jmQyBAvTcQCYN0wfzx8dB1ENYPkedDzGu/dX0REApJ69kRERM7C/mI3N7y5lG37i3ju2l7cNqitdxsozqP3+mchLwX63AxXvgJhUd5tQ0REApKKPRERkTOUsusQf112FBsUwgf3nMf5nZt458ZHD8Hu5bBrMWz4kkaF2TD6JRh4L5zNQuwiIlKvqNgTERE5jXK3ZUt2IVtzCknPKWRrThFbcwrZefAITSINkx8YRqdm0WfegLOEprmL4ZuZsGsJ5KQBFoLDoPW5rOn0KInn3ue15xH/evvttwHYsmULXbt29XMaEQlkKvZERESAl2dv5a15Rwia8+1Jx0qdbuysBQAEGWjXJIouCTGM6deSju59Z1foHc6EKbfSM2sthEZBm/Ng+FPQdii0SoTQCAodjjO/v9Q6ycnJADgcDpKSkvwbRkQCmoo9ERGp95wuNx8t3Umr6CBG9Gt30vGsfXu45NyedG4WQ4emUUSEBh875nBknnnDe1Ng6lgoKyKtx+/pef3vITj0zO8nIiJyHBV7IiJS7y1KP8ChYid3DQjn16O6n3Tc4cghqV8r7za6/nOY/jBEJ8CdX7F/Y44KPRER8SotvSAiIvXe9LX7iI0MpVeT4NOffLasG+b8Fb64t2KY5rh50OzkAlNERORsqWdPRETqtaNlLmZtzOHqfi0JCcrzbWOlRfRMewEOLIMBd8KolyAkzLdtiohIvaViT0RE6rU5m3MoLnNxVd+WlO3xQbHndlcuofAFbJxOk6P5cPnzMOgBLaNQT6WkpAAVs3HGxMSQmJjo50QiEqhU7ImISL329dpMmsWEM6h9PAv3eOmm1hJzOB2+mwVpX0JhVsVMm91GsSYkkQGDH/RSQ1IXDRw48IRta62fkohIoFOxJyIi9VbBUSeOLfu5fXBbgoO81MvmLIEPriRx78qKdfI6XQa9r4cul0NYFIe1jIKIiNQQFXsiIlJvfZ+WTZnLzZh+Lb130/nPw96VbOt4L51ueBoi47x3bxERkWrQbJwiIlJvfb02k7bxDejbOtY7N9y3Gha/Bv1vZ2+bMSr0RETEr1TsiYhIvZRbWMKSjANc1aclxhsTpZSXwfRHIKopjHju7O8nIiJyljSMU0RE6qVvUrNwW7w3hHPRy5CbBrdMUY+eiIjUCurZExGReunrdZl0ax5Dl4SYs79ZThos+Af0ugG6jTr7+4mIiHiBij0REal39he7Wb073yu9esbtgukPQ0QsXDHeC+lERES8Q8M4RUSk3lmeXQ7AVX3OvthrvXc6ZK6BGyZBVPxZ309ERMRb1LMnIiL1zrLMcgacE0ebxg3O7kYH0mm/YzJ0uxJ6XuudcCIiIl6iYk9EROqVrTmF7C2yjOl7lr16bhdMfwRXcBiMfgm8MaOniIiIF6nYExGReuXrtZkYYPTZDuGc/TTsWca2TvdBTHOvZBMREfEmFXsiIlJvWGuZkZpJj/ggmsaEn/mNVr4DS9+A85LJaX6x9wKKiIh4kYo9ERGpN7bmFLHzYDEDE85ifrL02fDN76DzSBj5d++FExER8TLNxikiIvXGrLRsAPo3Cz6zG2RvgGl3Q0JPuOE9CNavUam+cePGAZCZmUnLlmc/I6yIyKnot5SIiNQbszbm0P+cOOIinNW/+HAWTL4JwmPg1k8hPNr7AaVemDBhAgAOh4OkpCT/hhGRgKZhnCIiUi9k5h9l/b4CRvSo/mQqQa4SmHIzHM2HsZ9CbCsfJBQREfEuFXsiIlIv/LApB4DLeiRU70K3ix4bX4Ls9RVDN1v09UE6ERER71OxJyIi9cKstBw6NI2iU7NqDr+c/TRNDq6Ay5+Hrpf7JpyIiIgPqNgTEZGAV1DsZNn2g9UfwpnyASx9g72tRsOg+30TTkRExEdU7ImISMCbtyWXcrdlRM9qDOHcsQBm/gY6XkJGx3t9F05ERMRHNBuniIgEvNkbc2gaE06/1nGeXXAwAz69Axp3hBsnYZet8W1AqVeMMSdsW2v9lEREAp169kREJKCVOF04tuRyafcEgoLM6S84eqhiiQUTBGOnQkSs70OKiIj4gHr2REQkoC3NOMiRMpdHQziNuxw+uwsO7YK7vobGHWogoYiIiG+o2BMRkYA2a2M2UWHBDO0YX/WJ1tJp20TInA9X/xvaDq2ZgCIiIj6iYZwiIhKw3NYye2MuSV2bER4SXPXJKZNolfkdDHsU+t9WMwFFRER8SMWeiIgErO35bg4UlZ5+CGdpEcx9lkNxveGS/6uRbCIiIr6mYk9ERALW6lwXIUGGpK7Nqj5x1XtQfJAd7W+HIP1qFBGRwKDfaCIiErBW55QzpGM8sZGhpz6prBiWvAYdhnM4tlvNhRMREfExFXsiIhKQtuUWkV1sGdHjNEM4U96HI/vhosdrJJeIiEhN8VmxZ4x5zxiTa4zZcNy+xsaY2caY9MrvRpX7jTHmNWPMNmNMqjFmgK9yiYhI/TBrYzYAl1ZV7DlLYPGr0O4CaDukhpKJiIjUDF/27L0PXP6zfU8Ac6y1nYE5ldsAVwCdKz/JwJs+zCUiIvXArLQc2jcMokVs5KlPWvMRFGXDRb+vuWAiIiI1xGfFnrV2AZD3s91XAx9U/vwBcM1x+z+0FZYBccaYFr7KJiIigW1JxgHW7smnf0IVyy2Ul8KiV+CcIRU9eyIiIgHGWGt9d3Nj2gEzrLW9KrfzrbVxlT8b4JC1Ns4YMwN43lq7qPLYHOBxa+2qX7hnMhW9fyQkJCROnTrVZ/nPVFFREdHR0R7vP5NramMbgf58aqNutK026kbbvmxjc56Ll1NKaBJp+FVPN80b/fI1nQ8vouvWN1nX5/841Lh/rXuOmmwj0J+vtrUxfPjwE47PmzevxtpWG4HRRqA/X11rw9+GDx+eYq0d+IsHrbU++wDtgA3Hbef/7Pihyu8ZwPnH7Z8DDDzd/RMTE21tNG/evGrtP5NramMbgf58aqNutK026kbbvmpjxY6DtvufvrUX/2OezT1ccsprHHNmW/tyL2snXGyt233GbZ9t3trSRqA/X21rAzjhU5Ntq43AaCPQn6+uteFvwCp7inopxJdV5i/IMca0sNZmVQ7TzK3cvw9oc9x5rSv3iYiIeCRlVx53v7eC5g0jmDJuME1jwk95bkKOAwp2w+iXwJiaCykCDBhQMQ9dYWEhMTExfk4jIoGspou9r4G7gOcrv6cft/8RY8xUYBBQYK3NquFsIiJSR63ZfYi73ltJs4YRTEkeTLOGEac+2VVO213ToEU/6HxZzYUUqZSSkgKAw+EgKSnJv2FEJKD5rNgzxkwBkoAmxpi9wDNUFHmfGWPuBXYBN1We/g0wCtgGFAP3+CqXiIgElu0FLl55dwXx0WFMHjeIhKoKPYD104gsyYaLXlavnoiIBDSfFXvW2ltPceiSXzjXAg/7KouIiAQel9syIzWTf6wsIb5hJFPGDa56mQWA8jJY8CJFUe2J7jqqZoKKiIj4SU0P4xQRETkrbrflmw1Z/POHdLblFtEmJogp4wbTMu40hR7A8rcgL4Ptvf9EH/XqiYhIgFOxJyIidYLbWr5dX1HkbckppHOzaP41dgCRBzfTulGD09+gMAfmj4fOI8mL/+UZqkVERAKJij0REamVrLVkHy5hS3Yh6TlFfLikhD2Fq+nQNIpXb+nHlX1aEhxkcDi2eHbDOX+B8hK4/O+wfo9vw4uIiNQCKvZERKRWyDtSxqy0bL7dUMprGxeTnlNEYWn5sePNowyv3NyXMX1bERxUzSGY+1Jg7ccw9FcQ3xFQsSf+k5iYCPy09MKPs3OKiHibij0REfGb/OIyZqXl8N/UTJZkHMTltkSHQs/WQVzTvxVdEqLpnBBDl4QYUlcuIal/6+o34nbDt49DVDO48HfefwiRalq9erW/I4hIPaFiT0REatx3G7L5d0oJm2b/gNNladM4kuQLOzC6dwv2b13N8OFDvNfY+s9g70q45k2IaOi9+4qIiNRyKvZERKRGpew6xAMfpxAfYbhnWHuu7NOC3q1iMZWzYzrSvTdLZnB5Mcx+GlolQp9bvHZfERGRukDFnoiI1KhX56QTHxXG34aEMPLS7j5tq+2uaVCUA7dMhqAgn7YlIiJS2+g3n4iI1Jg1uw+xYOt+xl3YgfAQH69zdzCD1nu/hr5jobWWWhARkfpHxZ6IiNSY1+duo1GDUO4Y3Na3DR3Nh++ewB0UCpc+49u2REREaikN4xQRkRqRujefuZtz+d3IrkSFe/nXT9kR2L0UdiyAHQshay1YNzs7/g+dYpp7ty0REZE6QsWeiIjUiNfmbCM2MpQ7h3ihV684D/Yshz3L6Z/6HSzYBm4nBIVC63MrllhofxF7d5TR6exbExERqZNU7ImIiM9t2FfAD5ty+M1lXYiJCK3yXON2VhRzZUeO+xTRPGs2TJ8Gu5fDwfSKk4NCMFHtYegj0O4COGcwhEX9dLOdDt89lIiISC2nYk9ERHzu9bnpxESEcNfQdhWLnKfPovWe7+D72VCYBYXZx74vchbDgpPv0Q0gsjG0GQT9xlYUdi37s3rxcpKSkmr2gUREROoAFXsiIuJTm7IO831aDo9e0pnY7GUw6ynIWlcxvHJXJDRsATEtoGV/iG7O9pwCOnTtXdFDF9YAwqIhLIrlG3cx6IqxYHw8i6eIiEiAULEnIiI+9cbcbfQI388juR/D4m8gtg1cN5GFuVFccMmok4q33Q4HHQYnnXSfozvLVeiJiIhUg4o9ERHxmf2HCkjc9Cavhc4meFckXPI0DH4IQiNxORwq3kRERHxIxZ6IiHhf9gZYO5krUt8nMriY8t53EDziTxDdzN/JRERE6g0VeyIi4h1HDsL6abD2E8hOxQaFsrC8P5n9H+Xe6670dzqRWsNaC4DD4dDkQiLiUyr2RETk7GSvp+eGv8OClIq17lr0w335eO5edQ4pOW4WjLjU3wlFRETqJRV7IiJy5nYsgCljibUGBt1fsSRCQk8mLdrBgr0bSe4TTnx0uL9TioiI1Esq9kRE5MxsmgGf3wONO7Cq0+8ZOvJ6AHYeOMKL32/m4m7NGNKiyM8hRURE6q8gfwcQEZE6aM3H8Nkd0LwP3PMtZeHxALjdlse/SCU0KIi/Xdsbo9k2RURE/EbFnoiIVM/i12D6w9D+IrhzOjRofOzQJyt2s3xHHn+8sjvNYyP8GFJEREQ0jFNERDxjLR0yPoA9X0LPa+HatyHkp/fx9h4q5vlvNnFB5ybcNLCNH4OK1G7JyckAZGZmMnnyZCZMmODnRCISqFTsiYjIT6yFnYtokTkLFq2Fknw4egiO5kPBHs7ZlwKJ98DolyAo+LjLLE9+uR6Av1+n4ZsiVZk4ceIJ2yr2RMRXVOyJiEiFnI3w3eOwYwFdAbYCQaEQ2ajyE0dGh7voeOUr8LNibsG+chamH+CvV/ekdaMG/kgvIiIiP6NiT0SkvivOg3l/g1XvQnhDGPUPlubFMeTiURDa4ITCbo/DQcefFXp78oqZurmMQe0bc9ugtjWdXkRERE5BxZ6ISD1l3C5YMRHmPQclBTDwXhj+B2jQmFKHA8KiTnltfnEZ36dlMyM1iyUZBwk28ML1fQgK0vBNERGR2kLFnohIfZS5hsSU38CRndDuArjiBUjoWeUlhSVOFu518v6kFSxKP0C529I2vgEPXNSB1s5M2jU5dXEoIiIiNU/FnohIfeJ2weJ/wry/ERoSCzd9CN3HnPQO3s9t2FfAfR+sIvtwGa0bFXHvBe25sndLerVqiDEGhyO7hh5AREREPKViT0SkvsjfDV/eD7uXQM9rWRl3Pef3uOq0l81Ky+bRqWtp1CCUJ86L4P5rh2u2TRERkTpAi6qLiNQHqZ/Bm8Mgez1c8xbcMIny0JgqL7HW8vb8DO7/OIUuCdF89fAwujUOVqEnIiJSR6hnT0QkEJUWQVEOFGbTfeNLkLsA2vz/9s48TK6i3P+f07MlmZnsyYQEMHsgkJBNVoFJAogsoggKKu6g/lQU9Irr1auXe/WiV/SqV/EqbsgiIDsEkEzCThLISjKTfZ/JTDKT2Zfurt8fVcU509OdkHVmOt/P8/TT3W+fU9+qOm/VqbfqnNNnwpW/hUGj97t7PGm45YHl3Ld4G5dOOY6ffvA0+uTl8OaRz7kQQgghDhMK9oQQorfT3gQv/pxpSx+F5a3QUAUdTW/9PCzIgTnfgXNugpz9d/t1ze38ZHEra/Zs48Y54/nKBRP1lE0hhBCiF6JgTwgheivGwJrH4Klvwt6tBP1PgpHToajEvopHQFEJr66t4azzrk6bRFs8wYbqJiqqGtyrkTe21FHXlOT2D03jfdNHHeVCCSGEEOJwoWBPCCF6I7vXw5O3wLpnYPgp8MkneWNjO6WlpV02bdta1vl7PMG9i7by2xeaqXx6HomkASAnFjBmaCGnjxnEjH51CvSEEEKIXo6CPSGE6E10tDB6413w/EOQUwDv/k84/QZ7eebGsn3u2h5Pct/irfxq/jp27m1l/MAYnz9/LBNHFDOxpIgxQwspyM0BoKxs32kJIYQQouejYE8IIXoLFU/DE19jdN1mmHI1XPTv9lLN/RBPGu55bQv/89w6tte1MOPEgdx21Wl0bFvB7NmTjkLGhRBCCNEdKNgTQoieTv0OeOob8ObDMHQiS0/7IdPef+Pb2nXZ1jq++XwL1S0rOO2Egdz6/lM5f+Iw+0fo2/XQFSG6g8WLF7/1PmvWrG7OjRAim1GwJ4QQPZQgmYCXfw3zb4VkHOZ8F86+kboXXnpb+yeThm88uIJ4Ev7wiVnMnjRc/5EnRA9g5syZADQ0NLz1WQghjgQK9oQQoieybQkzXv8qNG6E8RfAJT+BwWMOKIlHlu1g9c56PndaAXNOKjlCGRVCCCFET0XBnhBC9CRaauGfP4DFd5KfPwiu/hNMvgIOcEWuLZ7gJ0+Xc8rI/pw+In6EMiuEEEKInoyCPSGE6AkYA8vvg6e/Dc274YzP8VreeZx7yiUHldxdr2xhW20L/3nlFBLbVx3mzAohhBCiNxDr7gwIIcQxT3UF/Oly+McNMPBEuKEM3vMjErn9Diq5htYOfjl/He8aP5RzJww7rFkVQgghRO9BK3tCCNFdNO+x/5m38B+Q3w8u/W+Y+QmI5RxSsr9buIE9Te3ccvFJhyefQojDyh133AFAeXk5FRUV3HDDDd2cIyFEtqJgTwghDgcdrdCyh5x4U+ZtjIHqcqh4CirmwdZXGG2SMPUauOiHUDT8kLOxq6GV3z2/kcumHseU4wcccnpCiMPPZz/72U7fFewJIY4UCvaEEOLtkkzA+vnw5j+YumkllAPNe+w9dh02yDsX4LViGHA8DBhl3/sfz/i1b8DSG6Fus01rxFQ492ssbjqOWZd/6rBl8Rf/XEtHIsnXLtKfpQshhBDHOgr2hBBif9RuhqV3wRt3Qf026DOQ3LzhUDQahp8M/YZA30HQdxDrVy9n3LA+sHebfe1YCs01HBfLh/Fz4V03wYSLbCAINJaVHbZsVjYlufu1rXzkjBMZPbTwsKUrhBBCiN6Jgj0hxLGNMeS118He7ZBog0QHJNoh0c7wqoXw59thQ5nddtwcePetMOkSXn/hJUpLS7skt7WpjHGp9o4WXnj+Bc6fc+ERLcoDa9spyI3xpTkTjqiOEEIIIXoHCvaEEMcu7U1w70c5Z/1z8FLXnycDDDgBSr8B0z4CA084OJ28vphY3qHkdL8s21rHosoEN86dwLDigiOqJYQQQojegYI9IcSxSeteuOuDsO01Nr3jQ4yeeg7k5ENOHuQWQE4+r6/ewIzLbzjkp2MeaVo7Enz34ZUU58P1547p7uwIIYQQooegYE8IcezRvAf+eiVUroCr7mRT9UBGzyztsln9jrIeH+glk4av/n0ZK7bv5YvTCijuc2RXEIUQQgjRe9Cfqgshji0ad8EfL4OqN+Gav8Ep7+vuHB0Stz9bwePLd3LLxScxs0Tzd0IIIYQIUbAnhDhmKGitgTsvgdqN8JH7YOK7uztLh8RDb2znF8+t44Ozjuez543t7uwIIYQQooehaWAhxLFB7SamLf0WJJvhow/CO87q7hwdEos37eHr9y/nzLGD+ff3TSEIgu7OkhBCCCF6GAr2hBDZSTIJlcth/XP2teUVcmMF8MlHYNTM7jJhOfwAACAASURBVM7dIbGrOcnNf1nCqEF9+c1HZ5Kfq4s0hBBCCNEVBXtCiOygowWqy6FyBSe/eR+89mlorrG/jZgCZ/0/Xu84iTN6eaBX39rB7a+3kkjm8vuPz2Jgv/zuzpIQQggheigK9oQQvY/WvQytfgnmvwS73oRdq2HPBjBJAAblDYSTL7Z/gj62FIpLAGgpK+u2LB8OmtrifP6vS6hqMvzlMzMYO6you7MkhBBCiB6Mgj0hRO8g3gZrn4bl90HFPE5NtEEQg8FjYfjJcOpV9r3kFF5asZXS2XO6O8eHlZ17W/j0HxezprKeT52az9njhnZ3loQQQgjRw1GwJ4TouSSTDKhbCY88CG8+ZP8Ivd9QmPkJXu8Yw4xLPgF5fbvuF2w/6lk9kizfVsdn/rSY5vYEv//4Owkq3+zuLAkhDoHf/va3AJSXlzNp0qRuzo0QIptRsCeE6FkYAzuXwcr7YeWDTK/fDnmFcPJlMOWD9rLMnFzqy8rSB3pZxpMrdnLTfUsZUljA/Z8/nZNG9KdMwZ4QvZobbrgBgLKyMkpLS7s3M0KIrEbBnhCiR9C3eTuU/QhW3A+710IsF8ZfwJujrmHy+78K+YXdncWjijGGx9a3c/9TrzP9xIHccd0shhUXdHe2hBBCCNGLULAnhOg+6rbCqn/Aygc4Y+dSIIDR74KzvgCTr4B+g9lVVsbkYyzQ60gk+cYDK3hgbQeXnzaS266aSp+8nO7OlhBCCCF6GQr2hBBHlfy2PfDKb2DlA7DtNWscOYN14z7J+Cu+Dv1Hdm8Gu5mORJIb736DJ1dWcsW4PG6/Zpr+MF0IIYQQB4WCPSHE4aWtgZLK+fDSSvtAlbdeddC4i7N2vAEYKJkCc/8VTnk/DB7LtrIyxh/jgV48keQr9yzlyZWVfPeyyYyLb1agJ4QQQoiDRsGeEOLwkOiA1/8EZT/i5KZqWAMQQJ8B4avvQDaNvoYxl94MwyZ2d457FPFEkpvuW8bjK3by7UtO5tPvGkNZ2ebuzpYQ4giwZMkSwD6Ns7i4mJkzZ3ZzjoQQ2YqCPSHEoWEMrHkcnv2+fbDKiWfzxoSbmX7xRyG/CGKxTptvLitjjAK9TiSShq/9fRmPLtvBLRefxPXnje3uLAkhjiCzZs3q9N0Y0005EUJkO7H9byKEEOnpv7cc7nwP3PsRCAK45m745BPsHTgZ+vTvEuiJriSN4V/uX8ZDS3fwL++exOdLx3V3loQQQgiRJWhlTwhxYLTUwpsPw/K/M2PzC1A4HC67HaZfBznqUg6EeCLJnSvbeX77dm6+cCJfmD2+u7MkhBBCiCxCIzMhxH6JJdph1UOw4u+w9mlItMOQ8WwYcx1jr/kRFBR1dxZ7BPFEkpa4oaG1o8tvVU1Jnl5VSUVVAxVVjVRUNbChuon2RJIb507gxrkTuiHHQgghhMhmFOwJITLTWA0Lb+PsJX+BRDMUlcA7PwNTroaR09myYAFjFegBsHxbHZ/9yxJ27m2FZ59Ov9Hz9qEMowb2ZWJJEedPGkZB/XZuukCBnhBCCCEOPwr2hBBdaW+GV34FL/wcOpqpGX4uIy76Mow5H2L6c+9Unlixk5vvW8qQwgI+NCmfCeO73ne3Y/N63nv+Oxk/vIiigrDrLSur0t8rCCGEEOKIoGBPCBFiEvDGX+G5W6FhB0y6FC74PmtW7WDEuNLuzl2PwxjDr8vWc9u8cmacOJA7PjaLlYtfpvTcrk/TLEtsYdoJA7shl0IIIYQ4VlGwJ4SAjhZY/xyzFn8TmjbDqJlw1e/hHWe7DXZ0a/Z6Ih1Jw9f+vpwHXt/GFdNG8uMPTKVPnlY9hRBCCNFzULAnxLGIScLOZbB+PmyYD1tegXgrOX1K4Ko/wClX2r9SEGnZ09TObYtaqajdxk0XTOTGueN1KaYQQgghehwK9oQ4FkgmYNdq2PoKbH6Js8ufhQV77W/DToZZn4Kxs3ltW8D5p17YvXntYRhjqKxvpaKqkbVVDVRUNfD82hqqG5L84trpvPe0kd2dRSGEEEKItCjYEyIbaW9mYO1yWLAItrwM2xZBW739raiE2kHTKDn7WhhbCv2Pe2s3s6OsGzLbMzDGUN3YxtqqRsorG1i7q4FFFS18qexpGlrjb203tCifk0b05zMno0BPCCGEED0aBXtCZAMttbDlVdj8og3udrzBtKQLUIZPhlM/ACeeCSecAYNGs3rBAkqmlXZrlruT3Y1tdqVul12pW1TewlcWPkNdc/j/eAP75TG8AK6YNpKJJcVvvQYX5gNQVlbWTbkXQgghhHh7KNgTorfS0QIv/4pZi/4CZZsBA7E8GDUDzv4Sy/cWMfWST0PfQd2d06NCRyLJi+tqeG7NLtZvbuPRXcs6/W4wrN7UwtdeeIaaxva37P375FLSB95z6nFMLCliYkkxE0qKGFZUwIIFCygtnXK0iyKEEEIIcVhQsCdEb8MYWPMYzPsW1G2hY+CpMPtb9smZo2ZCXl8A9pSVZX2gF08keWXDHh5bvoOnVlVS19xBYX4OfWJJNjXv7rJ9gYG5J5UwwQV1E0uKKemvoE4IIYQQ2YmCPSF6E9Xl8OQt9gmaw06Gjz3Csi2G0vNLuztnR4WmtjjrdjVSUdXA46va+Orz/2R3UzuF+TlcOLmES6eO5LyJQ3n5hecpLS3tsn9ZWRmlpVOPfsaFECLC9ddfD8COHTsYOVL3/gohjhwK9oToDbTuZdy6P8DCxyG/EN7zXzDr05CTC1vKujt3ACSShi17mqmoauC5jR2sZn3a7TZsaE/7Wyb7ivJ2/rxpERVVDWyrbXnLXpADF506kkunHEfppGH6jzshRK/hjjvuAPwEVGn3ZkYIkdUo2BOip7NjKdx9Lcc37ISZH4c534XCod2WnWTSsK22hYqqBp7a0M7DVUspr2xgfXUjbfFkuGH5msyJVGT4LY09N4DxJS1MP3EQH5p1AhNKiplYUsSmlYuYM3v6IZZGCCGEECJ7UbAnRE9mzePwwGeg3xBen3EbMy+//qhJG2PYXtfC2qpGntzYzqO7llFR1cC6XY20dCTe2m7kgN1MKCnmnPFDmFBSzKSSYraveYM5peelTXfhwoWcd17X3zLZX3phIXNmd7Vv0Z+YCyGEEELsEwV7QvREjIGXfwlPfxdGTodr76FhyeojJGWoqm+zK3UbO3iiZhkVVY2s29VIY1v4/3LDi6uZWFLMtaefyMSSIiaUFFNVsZRLLpzdJc3a9UHGyyrzc9L/lskeU1AnhBBCCHFQKNgToocRJOPw2E2w5E6YfAW87zeQ3w84tGDPGENNYztrqxp4ZlMH8x5cwdqqBsqrGjr9afiQwl1MLCnmAzNGuUsmi9m1dhmXXdQ1qCvbqEBMCCGEEKKnomBPiJ5ESx1TVvwQapfCu26COf8KsdgBJ7OnqZ01exJsfXkT5VUN9g/Eqxqo7fSn4TuZOLyY9542kkkjipkwvJjq9ct5b7qgbrOCOiGEEEKI3oaCPSF6CptegMduZmDdOnjvL2HGdfvdpaU9QUVtgu2vbmZtVSPllQ2s3dUQ+dPwVRT3yWViSTEXnzrirf+Wq16/nCsumk2Qcolk2VYFdUIIcaRJ7XuNMd2UEyFEtqNgT4juZucy+OcPYN2zUDyS5VO/z7T9BHot7Qn+8somfrNgA3ua2oGVFObnML6kmDknDWdiSTEtlRu46qJzGNG/T9egblusi00IIYQQQmQXCvaE6C52r4f5t8LKB6DPQLjwh3D69dS9+GrGXVo7Evz1lc38ZsEGahrbOHfCUGYUN3D1hWczckBfYrEwgCsr28JxA/oejZIIIYQQQogeiII9IQ4XxkBHC7TVQ2s9tNUzsHY5rI1Dot29OiDRxoSKR2Hhs5CTD+d+Dc7+EvQdmDHp1o4Ez2zu4OsvzmdXQxvnjB/Cby6YwazRgykrK+P4Qf2OYkGFEEIIIURvQMGeEFHi7QysXQEraqCpBpqq3asGmnczo64GVveFZLzT65ymOljYYr9HmAawrKvMcUEOzPoknPd1KC5hT1M7dzy5hrte3UxzW5zg6Sc6bZ80hqSB08cM5hfXTufMsUOOXB0IIYQQQoisQMGeEGBX3JbdDQtvY1rdljBAC2LQbygUDoN+g+nIGwADSyCWA7FciOVBLIdd1XWMGnsSFPSHPv2hYAD06c/SVRVMm3WGXcHzr9x8Xly8gnMvuJS65nb+b145d764keaOBJecehxBUw3veMeJnbIXENCvcRufv/JM3WsnhBBCCCHeFgr2xLFNIg7L74GFt0HtJhg5nZWjruXU86+0AV7fQZ3++mBFWRmlpaVdkllbVsaoNPa6HQVwwuld7PVmPf/9TAV3vrCRxvY4l045ji/PncCEkmLKysooLT2pyz5lZTsV6AkhhBBCiLeNgj3RMzEGkkkw/pUAkyQn3gJtjWl3Sf+bISfeDC11Lh3zVnollc/BL78CtRvhuNPg2nth4rupWbAAhncNtg6W1o4EG6qbqKhqcK9GXlzbTEt8LZdMGcGX505k0ojiw6YnhBBCCCEE9LBgLwiCi4GfAznA/xljftTNWRJHGmOgbitULoedy+175QpK926FBV03PxfghfRJZfotk/1kwIyYSs1lf2RF4VlUVDZRsXQZKze18MvVL6XV2Ls3/W+Z7Nuqm9k17ymS7i+UcmMBY4cVMmN4Lt++6iwmj+yfvjBCCCGEEEIcIj0m2AuCIAf4FXAhsA1YFATBI8aYN7s3Z8cIySTUb4Pqcti1GqrLOXVrBez8bZdNT62phu2/jqy62ddptbthQ/9wJS6ZeGtFbmZjI6wp6pyQgXN2b4IFDc4QwNAJcMIZbBp4FqPHjLf3zAWBe4+xfsNGxo0bR0cySU1jOzWNbdQ0tFHT2E5DYxP9+nV9KmVTcwt9C4tIEmAIMEEMQ8CbDf14bMc0mu5PAksAKOlfwMAcKMiLdUkHIC/Db5nso4pjfPCssUwsKWJiSTGjhxSSnxujrKxMgZ4QQgghhDii9JhgDzgdWGeM2QAQBME9wBVArwr2Xr33R+SvK2PJ4p91+S2/rS2tfV+/HS57pt8CDMe3bKVtwQ4KTOtb9r05g2k0xbQ07umSTjwep6lxL4YYJgjsOwFtHXGq25pJBjEMMZLkkgzyMcRoayugwOR3SWtX7HRqhk9ha8EEdhSMoT3WF9qheu8uhu0Y3mX77dUnUVdVwJY9zZ1Wy0YPLSQRNNM/t+vlkPVBA/1zutrb8xq4esqJTCgpYlJJMROGFzOgX567Z+7MtHWY6bd92yemTUsIIYQQQogjSWCM6e48ABAEwVXAxcaYz7jv1wFnGGO+mLLdDcANACUlJTPvueeeo57XfdG26A9MbFpEusdoGEhr39dvh8u+r9+qGMym4Hg2cjwbg1Fs4njqgyKSySSxWNfVqgO1H8w+meyBSTKyOJdRRbG3XiWFAbmxgMbGRoqKirrsc7jshzOtbNfI9vJli0a2l08avUP7WNSYPXt2p9/nz59/1LSlkR0a2V6+3qbR3cyePXuJMWZW2h+NMT3iBVyFvU/Pf78O+OW+9pk5c6bpicyfP/+A7AezT0/UyPbySaN3aEujd2hLo2dpZHv5epoGdv71rdfR1JZGdmhke/l6m0Z3Ayw2GeKl9Msw3cN24ITI9+OdTQghhBBCCCHEAdKT7tlbBEwIgmAMNsi7Bvhw92ZJCCGEEOLwMmPGDAAaGhooLtZf7wghjhw9JtgzxsSDIPgiMA/71wt/MMas6uZsCSGEEEIcVpYssU+Atg/xKu3ezAghspoeE+wBGGOeAJ7o7nwIIYQQQgghRG+nJ92zJ4QQQgghhBDiMKFgTwghhBBCCCGyEAV7QgghhBBCCJGFKNgTQgghhBBCiCykRz2gRQghhBAi25k5cyYQ/vWCfzqnEEIcbhTsCSGEEEIcRV5//fXuzoIQ4hhBl3EKIYQQQgghRBaiYE8IIYQQQgghshAFe0IIIYQQQgiRhSjYE0IIIYQQQogsRMGeEEIIIYQQQmQhCvaEEEIIIYQQIgtRsCeEEEIIIYQQWYiCPSGEEEIIIYTIQhTsCSGEEEIIIUQWomBPCCGEEEIIIbIQBXtCCCGEEEIIkYUo2BNCCCGEEEKILETBnhBCCCGEEEJkIYExprvzcNAEQVANbO7ufKRhKFBzAPaD2acnamR7+aTRO7Sl0Tu0pdGzNLK9fNmike3lk0bv0JZGz+MdxphhaX8xxuh1mF/A4gOxH8w+PVEj28snjd6hLY3eoS2NnqWR7eXLFo1sL580eoe2NNKn1VNfuoxTCCGEEEIIIbIQBXtCCCGEEEIIkYUo2Dsy3HGA9oPZpydqZHv5pNE7tKXRO7Sl0bM0sr182aKR7eWTRu/QlkYvolc/oEUIIYQQQgghRHq0sieEEEIIIYQQWYiCPSGEEEIIIYTIRrr7caBH4wWcAMwH3gRWAV929lOBPUAb0Ah8w9lvAeJAEkgAtzv7DUALYNxrD/Bl4CPA1oi9HbgPWAvURexJoNrpdUTsxul7eyJln9o0dgO0un3iKfbaDBrJDBo+z/E0dpNBoyPyStWIZ9BIROo0Vbs9g0amcsTT2A1Qn8HensHuNVK1zT7K4X9rTVO+RAZ7Mo1GIkO+4sCjadIxWH+qdOlF7auw/taeYm/G/h9MalpJZ2tNk1ZtBnsig0Yr9v8u02k0ZkirPVIvqe1gUxqN7U47Xd2uz1BX38tg3wpUpdHeCWxIo20yaCSwvtOeJq06wj4ktXxtaTT2AlsyaPjjl5rW3gzlW5PhONXvoxy702j4z+k0mt0xSecLe9Lsk8T6eSa/SldXcdL3C96e2m58+zuQ/iKeQdtk0E5GypJq93lITSsZeaXmtyONhvePffVjqWk1Z7D7tFI1fDqZ+sR0aWWqK3/eSHeczH403u5xyuQjrfs4Tpl8wWRIK5Pda6frv/elka4Okxk0mrHtM1WjJYO9nfC8lZpWuu0N0JDht0za0TaQyRdSt28j83ghnY8kyHz8ajJo1Lg8p9pfxvaJ6dJKVz5/LPblC6m2Wqed+lvDPsqR6ZhnGnc0Zchvawb7vnwhXZ68RkOG3xozHKeODBqZ6q8tUvZ0xzxduznQ/tuPD9O12XT1G9Xwn6N9c5s7vq0p2/vxysvAjsj3eETHf05g2+1SrL/UuO+twL2Et8/9LrL9o0B/Zy9w260DXgVG668X3j5x4KvGmMnAmcAXgiCYDHwe+KsxpgD4L+Drzv5O4PfGmBhwF/BFZz8fe8DeC/wJKAJuAgKgP/Bn9xvAB4CPA3nYxlOC7YiGAqXAU1hn+gZwP1AIfAJ4yKV3P9YBAqfzRcKTi983H/g7ocP8r7MPTNG4A3jWpTUvovGI0wDIAV4kdOwfR+rvfYQN9xm3TS5wXkRjG/CcS/cLacrheYywgd3p9s0FVjgNA/wEO1jOBR53Grj83e/SbQFecNv/xaUL0NflowPbaO/BHv884BXCjuwXbptc4HnCBuvtYAf/vkwrsI0Xt90PsY1yHbDM2WKu7LkuvRJXjpjL93Nu/zuBh519S6R8d2L9LQe4xO3rg8eXnMYAbGcTd/v8m6uDya4OAle329x7X+zxrHPb1wBPRLZrdZ93AQ+6dAdiB+y+s/uJe49FNMBOitS5eniA8E9GmyIaLS6tABtgeY08bLBg3D7/4eq9EHsS8hofd3Uw0tWz17gf69MxbKDi7d8BFrrPX3XlSrg0n3b2kU7bn6S+6fJR4vIYuH0+FannuyIaj2CDqRjwK6yfBNj277cpcvn1gymv0Q9Y4LY3zp7A9h8dkf1/4MoYc/v4OlwIvO7y1R/bpyRdPd/s0pxEZ18YiT1OxSl1dW0GjW86u2cdYYBzl3vvi53Q8hoLXV78iSpVI8D6oddYC/zT5de3g0ZXFz9wGjnYNuv7hf905chxek857bXYk2rM5XWJs2/CHnPvS5Uujx3Yfizu0toY0b6N0Cd/TxgsbSP0y5URjY2RMlUAi9z+67BtwLjfXnUavr/qcPndSdh/f4KwX3+KsG99jrB/W43tWw3wU0KfLsCeI/zkypedPR/rL21YH5nq9s1zafk6vAHrG2DPF17jQVcWzxciGlucLQfb1vxAyWvkYPtWX44/RsrRENFejB1E+TaxBzuobXLavm/d5myt2L4n6cr9q0j5rsX6SOC29cfpRkI/rHLpG2w//rdI+e53+7W7uvXH6X+cRiv2+O919jsi5bg+otGA9V3j6t9/Dlz9+fPQ/xK2pyZCf/upK0cfrG8/6TTWu+95rp5eJjz/LnBlKAZ+7fKQBL7vNIqw7dL72zVOo6979xprCX2hLqJRh/WZJPbYeg0/4daB9beX0mjkYtvKky6trVi/jhH6gG+zC109FWLb1F5XN5tcnga7NBuwx/ynzn4m1l+9L1zn9g1cuX35biLsHxdhJ+B9//IPQn5MGMSVOdsAbP/mfaHc7VuEPQd6jf8X0dji0jXAG1h/AHsMawn7Hu8L/eg8KXgv4Xm2KaJxG3aSNQ/rz/443Yftp5Put/WEbfbaiIZfLGgBZhCe61YQHqfnCc8RUV/4Pdb/vC/8n9NoiWjkY8cqfmz1b4TjId/+wZ43fbvZAyx39r3Y4+/PMZWEAdaLhP3Y/xL2rX8gbLMvuDImXDoNLt2d2HNM4Orvp+5zmztufbH92K+wC0Tt2AWdFcAst/+JbvvPAiOc3tex7XULcL8xZpqrm2ZgNnAlcAFwscvHOdj20YL1j39x9k8DtcaY8cDP6DwWP2iOiWDPGLPTGPO6+9yA7bBGYSv+P91mv8M6zihgCraDBDsIjTn7u4ClxphHgR8RHsjB2M7sG1gH841mlnt/CJiLdVCwHcNpWGcZiT0R4PJzhrPvAG539kasA8awjrfYFw240H2udNsk0mi0RjQmRzQ2RTTaXHli2IHJJYSDnmmRfAx0ZY5qJLGd5kvOfnqGcoBtMDHCgaHHRN7zI3V1jtMw2BULX44C4BSnUUvYUNqB6YSd32JXrjgw3GlXuf39cRrh3qtcuV9233e6tPwgZquzt2IDAz9Q/T5hgDjVvb9G52O+GTjJ1dUIbCcKttM9zX1uAP7dfQ5cvlvd93+470ngZMLB97ZIvX3M2ZoIV1cgHMzHscfjuYi9r7NXAt92Gj5ACZz9Mmy94jR8sH5cpHwfdulFT+RgO78ip1EV0QDrfwHh7Js/HuMIA7EBhP52EWGf9RyhL0x39jh2Ze7fnD0P2y5bgbuxxwRs5zrebV/vdDyTXR02AJ8jPDFfHdHowB5PX+48Z59I2DEb7MRKgB2AeI0A6wd+lnRSpHyDI+VL9XXvrxOArxAepz4uzQdcnn07+pDTqMdO/PjjNDlSjiCikRfR2IwdMHrK3barsW3THyevEcf2ey84+8UZNHIjGi3YiRx/IveBVwLrP74fm0oYrK2KlGOG+60aOyj0g9MGbL9Q7fLqg55GbLtbhW0f8yJ1ONDlrRGYQ+jr57r3ddjA2k8iLI5ojMBOuKVq98MGQD6A8b6YIFwZANs2oPPsMdgJx6nYvmcB4WDzJGz/WokNZHz52rB1m4MdOHmfAjjevT+AvZrFa8ScRiXwfqyf4rbxGruw5wNcml5jPfZqGc8o9/5gisZSp9Hi9vXlMBHtEcC3nL0FGIZtd7XYgNn3e8cRTpz5ySiwPu3LFxD6SBLbXiqxg33vh2DrvRLr0/5YGGCQ01iMPV7+N6+xHnt+8L7eP1KOpohGo9un0pXjRcK+zk9kNGL7K68xxL23RnR8Xr2v9yX09Qrsubna5flzkf1nuPw3EE5sgW3XPuCMRerKH49qrL/80dnXRTRysBO3MawveI0kdvLSa/g+KVVju9Nod2V6wtn3YttNJfa8+lln34Htp99wdfhIJK/9XX4ewPaVnlXu/QGsD/rj1Ep4nKoI/bAd2/4rsWMTHxDEsf1HDmEgiiuX11iP7RO8r/8tolEX0aggbE97CNuNX13zwU/UF4a69xpswOHrcEBE4w1s+wR7fLy/PYQdR8SwvubPK/cTnoshHLM8ANzq0gMbrEx1ZWgh9OmoRgk2EIlh+55i9/mBFI317n0j9vzoy9EKzHSfX4toDAbGYOu0EPiMs+/G9jFxrF/9gzBAm4HtXw22zfjxxBTCMXMR4STCauxCRQLbJ80mXJnz45tc7PhkBdZvP4Dtlww2SJ2NHe+e7urEYMeLYH0zNwiCiTgfMsa8gj12a7CLJ7jyPOQ+P+M0AK7ALiaBPWZzgyDw58CD5pgI9qIEQTAaOzh8FSgxxvgDVIB1rlT7dYQzs0MJnfcy7AE+AdsR5bp9Po1tgDFvxzb2UcBot+8tbr9m7GD9Pc5+UYp9rrMPwDbewL3OIpy9GIIdWOzGOps/pqkas519UIr9B87e19VL4OwDCTuDjzqNFmyDLXT2b7u0wDZ2P/C4NEM5Ym6bADvzdGWkTIWEs00DsCdpXPm8RnmkrvKwjc9rXBcpxzBsR7UQO9AFexz6RepqLrYhRn/rwHZk/qQ7NaJRjw3GwXYcH3ef+2Bng3wA80lXnmexx9wHcrNcORqwx+9WZw+c3a+olRFSSLhC+d1IXTURDkSuwHZOAXaQ7E++re4z2AFYv4jd16GflfP2P0c0oitcCfediIbBznb5Y7vc5d+vRH/P2XekaF9GOAD2WuuxwXqe22cd4UB8OuFxqnMavm68T+PscWydX+ls1c7uV2bOiuz7kitHO3aF3+flbmf3K1Fe72cRjf6EgUAc26bi2Nk+v7Kfg/UNsO3/IsIT4IOEA7ErCQf833MaSbe9P07RYH0QdiLF59fXWQXWF7zGQ5HyjSJsT3+LlOPThH7oZ3vjTruQMMh6073vxh4/71dPuM8xbHDifeGvdD4e1zp7MqIxD/jXSPmKsH1wFfBzwomgYmyb9oP90c4+ADsojWHbg2+zE1xa3n5dJJ1cV397sANKX4dFTiMHDPq4tgAAEqlJREFUO9iodPuMc/YBLt++fHMiGjHCCbdpKdqXOXuAnaApwPbRHyM8bn3dq4XOkxnDXfly3bbeF3IJ+8NREY2+2JXdAkJf8MFQofv9OOykhe+rCpxGnnsf6+ynEQaIgwh9uhg7e+0DvncR8vGIxqhIOea6tBqw7dX7W05Eo5ow2Cty9TUI2+79QMrXtde4gnAya0qK9uhIWoMidXW9+1xC6MOnY/sYnO50Qj88J1JX05xGCdbH5jj7+ZFyjMIGzV7Da7+I7Ssh7DsLsIPq7xH6ep7TKMROFFzk7IVOox92cOsDgeNdGb39KmcvcOXog+0zryBcac93GrnYfs/7rk+vP7bP8BPVIyIaCewKMO77DKdhsH7o62pSRGMS9viA7R+Px7b/YkJfaCFsg7uxE+lgz715wDucjg+4Y9hj5ccsH4yUYW7EPonQp+uddl9n9wHlWa6++mLHCr7cudhj5uvwnIjGJW77E12efLmvT9G4OlKH3r4UO56CsM/rg+0XPoHtmyD0hSHYvtT333Fs+/cavi8Yh/U3b78qkk6e09iOHX/5vmer2/4sbJsZEMnX8S5vZxMufIyPaEwA3u3sBdixkde4AhvAe3z/PZfQp2/G+nEr1ge/6ez+XOBXzv0qmF+t3+rqqNTZ/dgtH9vHXIEdfwfYcVM+tj/sQ+jrP8ceN5/HwS7dADu+9Vd/fNHpD8EGps9hff2dWB89xaVzndMehR0XVGP9+xrs1QFDgiDYir1C6meEE2OrXH7B+kq0X98KYIyJu/z489tBc0wFe0EQ+KX2rxhj6tPYm1Ps78E6Zap9NnbGIcCu0LSk2AcQXlrlmYRbkTHGXIB1xBj2BHg+thE/HLGfh3UWsJdV+BNVDvbEEx3Ut7vvZ7l902n4ge7dKRrHObu/zCHAdmaFhCfEnziNvdjZH7+Cs4FwBnkCttHh3tOVozpSjssJr5v2q1Tt2ODgUmyjM07Ta/i6MoRBkC/fV902W509B9tRnenK1Iq9fKMde0Ib5rY32EClA3uCbsEeJ7AnBK9hsAM4sIPBRe5zEbbx+6DJzwyOxc5wHw8YY8xsl1Yf7MpJAeGKRrN7XYId5Hl8+/SXhRKpK//7GYQnm53hrvTHdnR+titqL3WfqyP24kid+NVqsB3y8Mh3r9GB9VcfiPkTBc7Wz2lHy1OMvQwM7HH1M6KjsHXr24sfJNVjBxF+MORPhGDb05nucyLF7gck0RPOGMIV6jj2shOwJ75iwplbX85+WP8e7777VQ+wJxPfNioJT559CAeOW7GD+qi238fPNvfF+p4PqvygugN78pwb+e4H7jHgS+7zzoh9DLY+/PHwK6tg/XCi+/zZiH1jJE/bImm9h3DA7E82vqxNhP7mB1JgT1Q+uP1oxD4JOygEO6PtNS4mPN612GPrV9193wa23frL+C4jXNX2l/P1IxwQJbGTLImI/ebI9gmsPw9y2xusj9U6jX7ud+/L5YSrk/WE7WlVRCOXcCY6ulo/hHBWuo5wdXQ4tg/ylx35e2r2Yo+/72vqXPliWD8436XlL8cqwh4nr7ETewmzDz79hADYQX0ce0wHE66mvuTSGuzqx2vMi2iMJGzzVYSDsGGEQY6/xNtrXObLYYw5xaU10OU9J7KP12gl9HVfL8Ox5yC/ar/bveJYHzgP234M9hLNqLYPCm4hDCQuww6IjTsWScJJOz952Y69XDOGHaD7VReDbfPeR06NaPv7fLyGP8++FNG+NlLuOsKJ2hbCvgM63wuWJDx+xtVRM/Yc5f2zmnCluIRwAnIvdvXC1+MZhH27vw+u1ZXRt/88wgnCPMLgZk9Eo4AwaGzB+oy/DeMEwvZcHNH4MOEg1t8X7CcRfZ0Mct99G73c2f15bhf2vOX7bwjbTRPhxBvYttmB7X8/jF0JB3t1Uavb/sPYNgt2Isv3D35S26d/uyvfyYTndbB+1OHKexrhcTo7ReOdzr49Yv9YpNw1hPeK9cWea0a63/xlvh3YCQTffzfTWaPE2TcQHpNrCX0hjl3xNdgA0G8P1jc7CAM4r7HAafjnA/h+3Ze5Cdv2/bmukXAC9yqsv/nz/iLC81eMcAXue1i/2erqcAjheCgX29+C7b/9VT5gz1N9nAZYn/L9T7uz+zLWEF7aH70P7wPYPq8B298XYVfwwPZT9dg2NA67iuf7lnlum37YfuJVbBB4EaF/XIs9B4Ed9xYAFcaYE7Djz28Q8inCy0aLCa8oOTJ098NTjuJDWvKwB+vmiK0c2xHNwzpfecR+IdZBfhGxb8YuOa8nvAfjm9gZ8iS201+FdZwWZ08QdvDlzn4GnW98LnPvOwhvME5iZ/58kPFnwo6hhvBSOj9obqbzDacVdL752N+jkEmjmfAabBP53bhyJiMa/v691sh2ScLruf3lcek0ouXw91P5NH05fBl9wBw9Ae4mDPaiN2PXuPcWwgdX+MAotQ59p1JF5/L4fCUiGr4O/aWcxu3nLyPymv4hHf6ywNZImknsTFD0YRZrI79Fy/FmZJtk5HNVZPtExN4a2TcesZt92H1dpT5Mojqim0hJJ0nXPCUj+6S+KjKUY3fE1hax74zUT6qGP4aJFPuGDBprI7aovSFib4zYow8Bim6fru5StVN/r49sE83vrsjnRtKXL1WjkfTlaM+gEX34SSPpy9HxNjUyfd6Xv6VqpuYxnpKW70/8vUo+D9sIH4TVTGff836V+sAHf9L3D7Py9qUpx9j73E5CH2xL0fDb+Pbh231F5HglIxqvRbSj9VIR0d4U0fDl83WTrv9uiXxOEt73liAcCEbroxV7iZ3fPtqWH82g4bWjx8NEtk3V8Pey+XLsSdGLniOSkX29RrQcfqDu62oLoS/4/jPq06k+0hH5bUuKdrT+vY94vzLYc2D0OEXrwPud1/TH2l9i69tPbWSbqEZbRCP6oIv6yPa+jKm+Hk85JhUpdq/t+1Dv6x0p27fT2d/2RNL09kwa6XwhnUZqu6mOlKM9RcPvW5ei4X0hmaLxWORYpI4Jom3Ua0Qf0OF9wffD/nhU0LkOvbb3BW+PnqNbIului3z2vhDtJ3z7iGrsjWg0RuwbIts3RdJtimzTQld/S+cLjRk0vL2ZsF9I1fDHKTrmS6fh200yRaM5kmZ1iob/7PseX0e+3fg0/HH3970mUn5/hvAS96imn+j0CwT+mO2J6PvJiHqXxrOE/UNUIzru9fZ27OWWtZF6Wum2uRfbjhoI+/Za7EpwlauvWmw7rcdeupmDPRclgSUpsUkjNtB+zX2fB5zlPudi/T7QA1reBu56198Dq40x/x356RHs6sRq7EF72NnnY530Uawjefs8bCS/lXBF5hHC5d3t2IcODCNcqWvCzko9jp018DeL3o11gDXYwYLBzjg8RPjQgN1u+1rsLJOfPfEdTwf22t4m7OyAnz1uxQaU/mEf7YQrE1GNbRGNdux9En7W9z7CBjqD8L7BcsJrwR93GoFLZ5nTqCG8x2wj4YzFnpRyLCc8sa+OlMMHswns5X13uf3bCDvKeyLl2Et4I3+9s/sVDlz91WNvxDXO5i9L7MAeVz875Z9EabCzgr4O/f1PSTrf/+Vn5esIL/f014JXEvrJj91+YDsJf8N2RaQcFU7P39DsV+SSbh/fKd0fKceDkW1KCVcf/T13TVgfe8bZo0/z20A4C9WGPVH6gcQ/nd1gfcGvAp4X0fCdp8E+dOKBiMY8wmPrNVoJbzyPY2fkPUsIT2gLUjTA1un7IxrrsP4H9rJBXw9PEz5cojpi96vfCWybPofwBOHTbMMONLz2wy5PYGdd/XYLI9q/JLxvcj3hzfZVhLPgcewkkD8Jj45oLHHvLdj24zV2Rvb5j0g5/AxlDdYn/PZ+YO1XvKMa/sEojdgZ13QaZRGNXYT3OcUJL8P0/uof8vBIRMM/aKAJ2w96Da+bwM7Aew3fLvzT7dYR9gtPEA7W/uz2B1un/Z39T4Rt00+2GGw/ssDZ27CXWvkBxgLCmeU9WF+PEZ7EwfZJ/ph/jrDfa8bOPne48j9H2Lf6Fb+N2GPj77XNI2xr/iEAYAcHfiC4ibDf86u80X4vwPYt3k/muzrMwbZf72ON2MutfB36PDVjLy1P1Whzn33fuhp7iRHY1cL7ncZmOk+krSE8R/g+oZlwMix6HmomPA+B9Tdfjlqsj+QQDopqXdlXRPIP9lj5h2R4Wx3hoD96fvIDxiY6nwP9Az2SwH87zRyn6VeW2whXSiGc1GvDtk1/+Xy7e3n/vT+i4evmA4THqYPw3NJA55WQJwhXzpcSXklST+jrK7D9ay7hQBRsOysjvKTfn7cC7MM7PLmE/bpvTzHsce3vtil3GjiNe93ntWk0/P39/kEauPS9T7dENLzNYC+f8xpxwhXj1ohGNZ2vRHkFu5KTGjwuJ/SFLYTt5kan7VcX/QrR1yJ1WEP43IOH3Xdv91fsdGBXXiA8TjiNDVhfyCc8z/r+7L5IWj5/t2J9xPuOvxKljbBvM1h/85NWvm/19e4nw1+JlMNPJqRqtLo68/fvfyeNRgLbtr0vtkTy/PUUDb9qtSmi0UZ4v56fkPcafpwSx467mwhX3nw/0YptU7nYy5z98ajHnp+9RrFLz7cvf0XUHsJgbQWhj+USBt3POM1C7Plluvveir2KLRd7LD9KOEHoHyC13tXHdsKx8rlum3bC+1UT2LbSD3u58RqXn2ewYwP/ELp3Yq8+a8DekkMQBNGrnr4D/MZ9foRwZfYq4DnjIr9DortX3I7Sqt67CIOLpe51CeFDSNrcQVjhbNEZJz+zdAm2AURns9tdWk/R+VHAfsZuB53/qiE6Kxmd4fIvb0+m+S2dzexn+3Qa0XwcqkbqzPCR0Eidcdxf+gfzSl1xeDsameo3utrwdjQy1eFLhCtR0Vc7nWfOovnxs9/p9kn3qOrEPrZvSWM3hE/XStXemUHDnyxS7XUZtvcrrunyVZthn0zl2FcbSOdvSezgJF1a+9JozaCRbsUuTvjXKKm/NWQon798Jp12unT8MUz3W+sBaPj6yVTuPRl+23sAGqlaqbaODL9l6kNNBnum/iXT9mYf2gfa7+1r+8Olsa9Xpu2PRv+9P40DOU77KvOB7nOg9gMt3/7OHen8oC3Db2370MhUt5nKkKkc7RnS2tfxy7R9Op/2/eeBaGSLLxyoRqY6NGQeD+2vTtJtm2lMcjBjrkw+va/+LdO5PJM905Uj6fSjV0ikK0MD4S1Q0Xa0FTsZ5K9m81d2bcMGjD/Ejl3aCK9e24x9KFwtdrLil9gg8iRs7LGa9H+98E/CMUI98H1n74N96Nc67KLC2MMRB3lRIYQQQgghhBBZxDFxGacQQgghhBBCHGso2BNCCCGEEEKILETBnhBCCCGEEEJkIQr2hBBCCCGEECILUbAnhBBCCCGEEFmIgj0hhBBZQRAEQ4IgWOpelUEQbHefG4Mg+PUR0Pt+RGNlEATvPcD9fxAEwQUHsH1pEASP7X9LIYQQwpLb3RkQQgghDgfGmN3ANLCBGNBojPnJEZb9mTHmJ0EQnAw8HwTBcGNMcn87BUGQY4z51yOcNyGEEMc4WtkTQgiR1URXxNxq3J+CIHg+CILNQRBcGQTBfwVBsCIIgqeCIMhz280MgmBBEARLgiCYFwTBcfvSMMasxv7J7tAgCC4KguDlIAheD4Lg70EQFLk0NwVB8OMgCF4Hrg6C4I9BEFzlfpsbBMEbLh9/CIKgwNkvDoJgjdvnyiNXS0IIIbIRBXtCCCGONcYBc4D3An8F5htjpgAtwKUu4Psf4CpjzEzgD8Ct+0owCIIzgCRggO8AFxhjZgCLgZsjm+42xswwxtwT2bcP8EfgQy4fucDnnf13wOXATGDEoRZcCCHEsYUu4xRCCHGs8aQxpiMIghVADvCUs68ARgOTgFOBZ4IgwG2zM0NaNwVB8FGgAfgQcAYwGXjR7ZsPvBzZ/t40aUwCNhpjKtz3PwFfAMqcfS1AEAR/BW44wLIKIYQ4hlGwJ4QQ4lijDcAYkwyCoMMYY5w9iT0vBsAqY8xZbyOtn0XvCwyC4HLgGWPMtRm2bzqEfAshhBAHhC7jFEIIITpTDgwLguAsgCAI8oIgOOVt7vsKcE4QBOPdvoVBEEx8G3qj/T7AdcACYI2zj3P2TAGkEEIIkRYFe0IIIUQEY0w7cBXw4yAIlgFLgbPf5r7VwCeAu4MgWI69hPOk/ezTCnwS+Lu7tDQJ/MbZbwAedw9o2XVwJRJCCHGsEoRXrwghhBBCCCGEyBa0sieEEEIIIYQQWYiCPSGEEEIIIYTIQhTsCSGEEEIIIUQWomBPCCGEEEIIIbIQBXtCCCGEEEIIkYUo2BNCCCGEEEKILETBnhBCCCGEEEJkIf8fP3Mh9oR4DrgAAAAASUVORK5CYII=\n",
            "text/plain": [
              "<Figure size 1080x720 with 1 Axes>"
            ]
          },
          "metadata": {
            "tags": [],
            "needs_background": "light"
          }
        }
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {
        "id": "EXEZPdSfeXXC",
        "colab_type": "text"
      },
      "source": [
        "## % change "
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "QTpEv8gEeXXD",
        "colab_type": "code",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 232
        },
        "outputId": "d9fabe97-40f2-43fb-945d-7a063103e0a5"
      },
      "source": [
        "data={'date':testYears, 'predictions':pred}\n",
        "pred=pd.DataFrame(data)\n",
        "act=pd.DataFrame(actual).reset_index()\n",
        "tot=pred.merge(act, on='date')"
      ],
      "execution_count": 202,
      "outputs": [
        {
          "output_type": "error",
          "ename": "NameError",
          "evalue": "ignored",
          "traceback": [
            "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m",
            "\u001b[0;31mNameError\u001b[0m                                 Traceback (most recent call last)",
            "\u001b[0;32m<ipython-input-202-4b258a81ec38>\u001b[0m in \u001b[0;36m<module>\u001b[0;34m()\u001b[0m\n\u001b[1;32m      1\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m      2\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m----> 3\u001b[0;31m \u001b[0mdata\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0;34m{\u001b[0m\u001b[0;34m'date'\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0mtestYears\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;34m'predictions'\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0mpred\u001b[0m\u001b[0;34m}\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m      4\u001b[0m \u001b[0mpred\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0mpd\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mDataFrame\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mdata\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m      5\u001b[0m \u001b[0mact\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0mpd\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mDataFrame\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mactual\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mreset_index\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
            "\u001b[0;31mNameError\u001b[0m: name 'pred' is not defined"
          ]
        }
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "Vu_m2IEPxV4l",
        "colab_type": "code",
        "colab": {
          "base_uri": "https://localhost:8080/",
          "height": 238
        },
        "outputId": "f30883eb-04a4-4cc5-f288-7495a915f314"
      },
      "source": [
        "actual"
      ],
      "execution_count": 282,
      "outputs": [
        {
          "output_type": "execute_result",
          "data": {
            "text/plain": [
              "date\n",
              "2020-03-10      0.0\n",
              "2020-03-11      0.0\n",
              "2020-03-12      0.0\n",
              "2020-03-13      0.0\n",
              "2020-03-14      0.0\n",
              "              ...  \n",
              "2020-07-06    280.0\n",
              "2020-07-07    281.0\n",
              "2020-07-08    285.0\n",
              "2020-07-09    294.0\n",
              "2020-07-10    299.0\n",
              "Name: Marion, Length: 123, dtype: float64"
            ]
          },
          "metadata": {
            "tags": []
          },
          "execution_count": 282
        }
      ]
    },
    {
      "cell_type": "code",
      "metadata": {
        "id": "QADdQ8G9eXXG",
        "colab_type": "code",
        "colab": {}
      },
      "source": [
        "change = round(100*\n",
        "      (\n",
        "          (pred.iloc[-1][\"predictions\"]/\n",
        "            act.iloc[-1][county_of_interest])\n",
        "          -1), \n",
        "      2)\n",
        "\n",
        "if change < 0:\n",
        "  word = \"fewer\"\n",
        "else:\n",
        "  word = \"more\"\n",
        "\n",
        "if good_policy == True:\n",
        "  phrase = \"a more lenient\"\n",
        "else:\n",
        "  phrase = \"a stricter\"\n",
        "\n",
        "\n",
        "print((\"As of {} (the date is in \"+format+\" form), {} County would have had {}% \"+\n",
        "       word+\" cases if it had followed \"+phrase+\" policy toward Covid-19.\").\n",
        "      format(\n",
        "          act.iloc[-1][\"date\"],\n",
        "          county_of_interest,\n",
        "          change\n",
        "          ))"
      ],
      "execution_count": null,
      "outputs": []
    }
  ]
}